Scippy

SCIP

Solving Constraint Integer Programs

reduce_base.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 reduce_base.c
17  * @brief Reduction tests for Steiner problems
18  * @author Daniel Rehfeldt
19  *
20  * This file includes several packages of reduction techniques for different Steiner problem variants.
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 /*lint -esym(750,REDUCE_C) -esym(766,stdlib.h) -esym(766,string.h) */
28 #define REDUCE_C
29 #define STP_RED_EXTENSIVE FALSE /**< just for testing */
30 #define STP_REDBOUND_SDSP 200 /**< visited edges bound for SDSP test */
31 #define STP_REDBOUND_SDSP2 800 /**< visited edges bound for SDSP test */
32 #define STP_REDBOUND_PCBDK 80 /**< visited edges bound for PC BDK test */
33 #define STP_REDBOUND_PCBDK2 400 /**< visited edges bound for PC BDK test */
34 
35 #define STP_REDBOUND_BDK 400 /**< visited edges bound for BDK test */
36 #define STP_REDBOUND_BDK2 600 /**< visited edges bound for BDK test */
37 #define STP_REDBOUND_SDSTAR 400 /**< visited edges bound for SD Star test */
38 #define STP_REDBOUND_SDSTAR2 800 /**< visited edges bound for SD Star test */
39 #define STP_RED_MWTERMBOUND 400
40 #define STP_RED_EXPENSIVEFACTOR 2
41 #define STP_RED_GLBFACTOR 3
42 #define STP_RED_EDGELIMIT 100000
43 #define STP_RED_EDGELIMIT_HUGE 1000000
44 #define STP_RED_MAXNINNROUNDS 15
45 #define STP_RED_MAXNOUTROUNDS_PC 4
46 #define STP_RED_MAXNOUTROUNDS_MW 4
47 
48 #define STP_BND_THRESHOLD 0.03
49 #define USE_FULLSEPA
50 
51 #include "mincut.h"
52 #include <stdio.h>
53 #include <stdlib.h>
54 #include <string.h>
55 #include <assert.h>
56 #include "graph.h"
57 #include "reduce.h"
58 #include "heur_tm.h"
59 #include "misc_stp.h"
60 #include "scip/scip.h"
61 #include "probdata_stp.h"
62 #include "prop_stp.h"
63 #include "portab.h"
64 #include "dpterms.h"
65 #include "relax_stpenum.h"
66 
67 
68 
71 
72 
73 /** returns limit parameter for SPG method */
74 static
76  const GRAPH* g,
77  int roundnumber,
78  SCIP_Bool fullreduce,
79  enum STP_REDTYPE redtype
80 )
81 {
82  int limit = -1;
83 
84  assert(roundnumber >= 0);
85  switch (redtype)
86  {
87  case stp_bdk:
88  limit = (roundnumber == 0) ? STP_REDBOUND_BDK : STP_REDBOUND_BDK2;
89  break;
90  case stp_sdstarbot:
91  limit = (roundnumber == 0) ? STP_REDBOUND_SDSTAR : STP_REDBOUND_SDSTAR2;
92  break;
93  case stp_sdstar:
94  limit = STP_REDBOUND_SDSTAR;
95  break;
96  default:
97  assert(0);
98  }
99 
100  if( fullreduce )
101  {
102  if( stp_bdk == redtype )
103  limit /= 2;
104  }
105  else
106  {
107  limit /= 3;
108  }
109 
110  assert(limit >= 0);
111 
112  return limit;
113 }
114 
115 
116 /** returns limit parameters for PCSTP method */
117 static
119  const GRAPH* g,
120  int roundnumber,
121  SCIP_Bool beFast,
122  enum PC_REDTYPE redtype
123 )
124 {
125  const int nedges = g->edges;
126  const int sqrtnedges = (int) sqrt(nedges);
127  int limit = 0;
128 
129  assert(roundnumber >= 0);
130  assert(nedges >= 0);
131  assert(sqrtnedges >= 0 && sqrtnedges <= nedges);
132 
133  // todo try mabe 100 500?
134  // or somethin similar
135 
136  switch (redtype)
137  {
138  case pc_sdc:
139  limit = (roundnumber > 0) ? STP_REDBOUND_SDSP2 : STP_REDBOUND_SDSP;
140  break;
141  case pc_sdstar:
142  limit = (roundnumber > 0) ? STP_REDBOUND_SDSP2 : STP_REDBOUND_SDSP;
143  break;
144  case pc_sdstar2:
145  limit = (roundnumber > 0) ? STP_REDBOUND_SDSP : STP_REDBOUND_SDSP / 2;
146  break;
147  case pc_sdw1:
148  limit = (roundnumber > 0) ? STP_REDBOUND_SDSP2 : STP_REDBOUND_SDSP / 2;
149  break;
150  case pc_sdw2:
151  limit = (roundnumber > 0) ? STP_REDBOUND_SDSP : 0;
152  break;
153  case pc_bd3:
154  limit = (roundnumber > 0) ? STP_REDBOUND_PCBDK2 : STP_REDBOUND_PCBDK;
155  break;
156  default:
157  assert(0);
158  }
159 
160  // todo more or less random number, tune them at least a bit
161  if( beFast )
162  {
163  limit /= 2;
164  }
165  else
166  {
167  limit = (int) MAX(limit, limit * sqrtnedges / 200.0);
168  }
169 
170  assert(limit >= 0);
171 
172  return limit;
173 }
174 
175 
176 /** print reduction information */
177 static
180  const char* method,
181  int nelims
182 )
183 {
184  assert(nelims >= 0);
185 
186 #ifdef STP_PRINT_STATS
187  if( print )
188  printf("%s: %d \n", method, nelims);
189 #endif
190 
191 }
192 
193 
194 /** iterate NV and SL test while at least minelims many contractions are being performed */
195 static
197  SCIP* scip, /**< SCIP data structure */
198  const int* edgestate, /**< for propagation or NULL */
199  GRAPH* g, /**< graph */
200  PATH* vnoi, /**< Voronoi */
201  SCIP_Real* nodearrreal, /**< array */
202  SCIP_Real* fixed, /**< offset */
203  int* edgearrint, /**< array */
204  int* vbase, /**< array */
205  int* neighb, /**< array */
206  int* distnode, /**< array */
207  int* solnode, /**< array */
208  STP_Bool* visited, /**< array */
209  int* nelims, /**< number of eliminations */
210  int minelims /**< minimum number of eliminations */
211  )
212 {
213  int elims;
214  int nvelims;
215  int slelims;
216  int degelims;
217  int totalelims;
218 
219  assert(g != NULL);
220  assert(vbase != NULL);
221  assert(vnoi != NULL);
222  assert(nodearrreal != NULL);
223  assert(visited != NULL);
224  assert(minelims >= 0);
225 
226  *nelims = 0;
227  totalelims = 0;
228 
229  do
230  {
231  elims = 0;
232  degelims = 0;
233 
234  /* NV-reduction */
235  SCIP_CALL( reduce_nvAdv(scip, edgestate, g, vnoi, nodearrreal, fixed, edgearrint, vbase, distnode, solnode, &nvelims) );
236  elims += nvelims;
237 
238  SCIPdebugMessage("NV-reduction (in NVSL): %d \n", nvelims);
239 
240  /* SL-reduction */
241  SCIP_CALL( reduce_sl(scip, edgestate, g, vnoi, fixed, vbase, neighb, visited, solnode, &slelims) );
242  elims += slelims;
243 
244  SCIPdebugMessage("SL-reduction (in NVSL): %d \n", slelims);
245 
246  /* trivial reductions */
247  if( elims > 0 )
248  {
249  if( g->stp_type == STP_PCSPG || g->stp_type == STP_RPCSPG )
250  SCIP_CALL( reduce_simple_pc(scip, edgestate, g, fixed, &degelims, NULL, solnode) );
251  else
252  SCIP_CALL( reduce_simple(scip, g, fixed, solnode, &degelims, NULL) );
253  }
254  else
255  {
256  degelims = 0;
257  }
258 
259  elims += degelims;
260 
261  SCIPdebugMessage("Degree Test-reduction (in NVSL): %d \n", degelims);
262 
263  totalelims += elims;
264  }while( elims > minelims );
265 
266  *nelims = totalelims;
267 
268  assert(graph_valid(scip, g));
269 
270  return SCIP_OKAY;
271 }
272 
273 static
275  SCIP* scip, /**< SCIP data structure */
276  GRAPH* g, /**< graph data structure */
277  PATH* vnoi, /**< Voronoi data structure */
278  int* heap, /**< heap array */
279  int* state, /**< array to store state of a node during Voronoi computation*/
280  int* vbase, /**< Voronoi base to each node */
281  int* nodesid, /**< array */
282  int* nodesorg, /**< array */
283  int* nelims, /**< pointer to store number of eliminated edges */
284  int redbound, /**< reduction bound */
285  SCIP_Bool verbose, /**< be verbose? */
286  SCIP_Bool* rerun /**< use again? */
287 )
288 {
289  // todo test properly
290 #ifdef SCIP_DISABLED
291  SCIP_CALL( reduce_sdPc(scip, g, vnoi, heap, state, vbase, nodesid, nodesorg, nelims) );
292 
293  if( verbose )
294  printf("pc_SD eliminations: %d \n", *nelims);
295 
296  if( *nelims <= redbound )
297  *rerun = FALSE;
298 #else
299  *rerun = FALSE;
300 #endif
301 
302  return SCIP_OKAY;
303 }
304 
305 
306 static
308  SCIP* scip, /**< SCIP data structure */
309  GRAPH* g, /**< graph structure */
310  PATH* pathtail, /**< array for internal use */
311  PATH* pathhead, /**< array for internal use */
312  int* heap, /**< array for internal use */
313  int* statetail, /**< array for internal use */
314  int* statehead, /**< array for internal use */
315  int* memlbltail, /**< array for internal use */
316  int* memlblhead, /**< array for internal use */
317  int* nelims, /**< point to return number of eliminations */
318  int limit, /**< limit for edges to consider for each vertex */
319  SCIP_Real* offset, /**< offset */
320  int redbound, /**< reduction bound */
321  SCIP_Bool verbose, /**< be verbose? */
322  SCIP_Bool* rerun /**< use again? */
323 )
324 {
325  SCIP_CALL( reduce_bd34(scip, g, pathtail, pathhead, heap, statetail, statehead, memlbltail,
326  memlblhead, nelims, limit, offset) );
327 
328  if( verbose )
329  printf("pc_BDk eliminations: %d \n", *nelims);
330 
331  if( *nelims <= redbound )
332  *rerun = FALSE;
333 
334  return SCIP_OKAY;
335 }
336 
337 static
339  SCIP* scip, /**< SCIP data structure */
340  SCIP_Bool usestrongreds, /**< allow strong reductions? */
341  GRAPH* g, /**< graph */
342  PATH* vnoi, /**< Voronoi data structure */
343  SCIP_Real* nodearrreal, /**< array */
344  SCIP_Real* fixed, /**< offset */
345  int* edgearrint, /**< array */
346  int* vbase, /**< array */
347  int* neighb, /**< array */
348  int* distnode, /**< array */
349  int* solnode, /**< array */
350  STP_Bool* visited, /**< array */
351  int* nelims, /**< number of eliminations */
352  int redbound, /**< reduction bound */
353  SCIP_Bool verbose, /**< be verbose? */
354  SCIP_Bool* rerun /**< use again? */
355 )
356 {
357  // todo propagate usestrongreds and use it!
358  SCIP_CALL( execNvSl(scip, NULL, g, vnoi, nodearrreal, fixed, edgearrint, vbase, neighb,
359  distnode, solnode, visited, nelims, redbound) );
360 
361  if( verbose )
362  printf("pc_NVSL eliminations: %d \n", *nelims);
363 
364  if( *nelims <= redbound / 2 )
365  *rerun = FALSE;
366 
367  return SCIP_OKAY;
368 }
369 
370 static
372  SCIP* scip, /**< SCIP data structure */
373  GRAPH* graph, /**< graph data structure */
374  PATH* vnoi, /**< Voronoi data structure */
375  SCIP_Real* radius, /**< radius array */
376  SCIP_Real* offset, /**< pointer to the offset */
377  int* heap, /**< heap array */
378  int* state, /**< array to store state of a node during Voronoi computation*/
379  int* vbase, /**< Voronoi base to each node */
380  int* nelims, /**< pointer to store number of eliminated edges */
381  int redbound, /**< reduction bound */
382  SCIP_Bool verbose, /**< be verbose? */
383  SCIP_Bool* rerun /**< use again? */
384 )
385 {
386  SCIP_Real ub = -1.0;
387 
388  SCIP_CALL( reduce_bound(scip, graph, vnoi, radius, offset, &ub, heap, state, vbase, nelims) );
389 
390  if( verbose )
391  printf("pc_BND eliminations: %d \n", *nelims);
392 
393  if( *nelims <= redbound )
394  *rerun = FALSE;
395 
396  return SCIP_OKAY;
397 }
398 
399 
400 
401 /** returns pointer */
402 static
404  REDBASE* redbase /**< base */
405  )
406 {
407  assert(redbase);
408  assert(redbase->redsol);
409 
410  return reduce_solGetOffsetPointer(redbase->redsol);
411 }
412 
413 
414 /** returns pointer
415  * todo: remove once redsol is properly tested */
416 static
418  REDSOLLOCAL* redsollocal, /**< data structure for retaining primal solution */
419  REDBASE* redbase /**< base */
420  )
421 {
422  REDSOL* redsol;
423  assert(redbase && redsollocal);
424 
425  redsol = redbase->redsol;
426  assert(redsol);
427 
428  if( reduce_solUsesNodesol(redsol) )
429  {
430  assert(!redbase->solnode);
431  assert(reduce_sollocalUsesNodesol(redsollocal));
432 
433  return reduce_sollocalGetSolnode(redsollocal);
434  }
435 
436  return redbase->solnode;
437 }
438 
439 /** inner STP reduction loop */
440 static
442  SCIP* scip, /**< SCIP data structure */
443  SCIP_RANDNUMGEN* randnumgen, /**< generator */
444  GRAPH* g, /**< graph data structure */
445  REDSOLLOCAL* redsollocal, /**< data structure for retaining primal solution */
446  REDBASE* redbase, /**< parameters */
447  SCIP_Bool* wasDecomposed /**< pointer to mark whether to exit early */
448 )
449 {
450  const RPARAMS* const redparameters = redbase->redparameters;
451  BIDECPARAMS* const bidecompparams = redbase->bidecompparams;
452  SCIP_Real timelimit;
453  SCIP_Bool rerun = TRUE;
454  int inner_rounds = 0;
455  int inner_restarts = 0;
456  const int reductbound = redparameters->reductbound;
457  const SCIP_Bool fullreduce = redparameters->fullreduce;
458  const SCIP_Bool nodereplacing = redparameters->nodereplacing;
459  const SCIP_Bool extensive = STP_RED_EXTENSIVE;
460  const SCIP_Bool usestrongreds = redparameters->usestrongreds;
461  SCIP_Bool sd = TRUE;
462  SCIP_Bool sdbiased = redparameters->usestrongreds;
463  SCIP_Bool sdc = FALSE;
464  SCIP_Bool pathrep = usestrongreds && redparameters->dualascent;
465  SCIP_Bool sdstar = TRUE;
466  SCIP_Bool da = redparameters->dualascent;
467  SCIP_Bool bdk = redparameters->nodereplacing;
468  SCIP_Bool bred = redparameters->boundreduce;
469  SCIP_Bool nvsl = redparameters->nodereplacing;
470 
471  *wasDecomposed = FALSE;
472 
473  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
474 
475  /* inner reduction loop */
476  while( rerun && !SCIPisStopped(scip) )
477  {
478  SCIP_Bool skiptests = FALSE;
479  int danelims = 0;
480  int sdnelims = 0;
481  int sdcnelims = 0;
482  int sdstarnelims = 0;
483  int bdknelims = 0;
484  int nvslnelims = 0;
485  int brednelims = 0;
486  int degtnelims = 0;
487  int sdbiasnelims = 0;
488  int pathrepnelims = 0;
489 
490  if( bidecompparams && bidecompparams->newLevelStarted )
491  {
492  skiptests = TRUE;
493  bidecompparams->newLevelStarted = FALSE;
494  }
495 
496  assert(*wasDecomposed == FALSE);
497 
498  if( SCIPgetTotalTime(scip) > timelimit )
499  break;
500 
501  if( sd && !skiptests )
502  {
503  SCIP_CALL( reduce_sd(scip, g, redbase, &sdnelims) );
504 
505  if( sdnelims <= reductbound && !extensive )
506  sd = FALSE;
507 
508  reduceStatsPrint(fullreduce, "sd", sdnelims);
509 
510  if( SCIPgetTotalTime(scip) > timelimit )
511  break;
512  }
513 
514  if( sd && !skiptests )
515  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &degtnelims, NULL));
516 
517  if( sdstar && !skiptests )
518  {
519  SCIP_CALL( reduce_sdStarBiased(scip, getWorkLimitsStp(g, inner_rounds, fullreduce, stp_sdstar),
520  usestrongreds, g, &sdstarnelims) );
521 
522  if( sdstarnelims <= reductbound && !extensive )
523  sdstar = FALSE;
524 
525  reduceStatsPrint(fullreduce, "sdstar", sdstarnelims);
526  }
527 
528  if( sdc && !skiptests )
529  {
530  SCIP_CALL( reduce_sdsp(scip, g, redbase->vnoi,
531  redbase->heap, redbase->state, redbase->vbase, redbase->nodearrint,
532  redbase->nodearrint2, &sdcnelims,
533  ((inner_rounds > 0) ? STP_REDBOUND_SDSP2 : STP_REDBOUND_SDSP), usestrongreds));
534 
535  if( sdcnelims <= reductbound && !extensive )
536  sdc = FALSE;
537 
538  reduceStatsPrint(fullreduce, "sdsp", sdcnelims);
539 
540  if( SCIPgetTotalTime(scip) > timelimit )
541  break;
542  }
543 
544  if( (sdc || sdstar) && !skiptests )
545  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &degtnelims, NULL));
546 
547  if( pathrep )
548  {
549  SCIP_CALL( reduce_pathreplace(scip, g, &pathrepnelims) );
550  if( pathrepnelims <= reductbound )
551  pathrep = FALSE;
552 //printf("pathrep=%d \n", pathrepnelims);
553 
554  reduceStatsPrint(fullreduce, "pathrep", pathrepnelims);
555  }
556 
557  if( bdk && !skiptests )
558  {
559  SCIP_CALL( reduce_bdk(scip, getWorkLimitsStp(g, inner_rounds, fullreduce, stp_bdk), g, &bdknelims) );
560 
561  if( bdknelims <= STP_RED_EXPENSIVEFACTOR * reductbound && !extensive )
562  bdk = FALSE;
563  else
564  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &degtnelims, NULL));
565 
566  reduceStatsPrint(fullreduce, "bdk", bdknelims);
567 
568  if( SCIPgetTotalTime(scip) > timelimit )
569  break;
570  }
571 
572  if( sdbiased && !skiptests )
573  {
574  SCIP_CALL( reduce_impliedProfitBased(scip, getWorkLimitsStp(g, inner_rounds, fullreduce, stp_sdstarbot), g,
575  redbaseGetSolnode(redsollocal, redbase), redbaseGetOffsetPointer(redbase), &sdbiasnelims) );
576 
577  if( sdbiasnelims <= reductbound && !extensive )
578  sdbiased = FALSE;
579 
580  reduceStatsPrint(fullreduce, "sdbiasnelims", sdbiasnelims);
581 
582  if( SCIPgetTotalTime(scip) > timelimit )
583  break;
584  }
585 
586  if( nvsl )
587  {
588  SCIP_CALL( execNvSl(scip, NULL, g, redbase->vnoi,
589  redbase->nodearrreal, redbaseGetOffsetPointer(redbase), redbase->edgearrint,
590  redbase->vbase, redbase->nodearrint, NULL,
591  redbaseGetSolnode(redsollocal, redbase),
592  redbase->nodearrchar, &nvslnelims, reductbound));
593 
594  if( nvslnelims <= reductbound && !extensive )
595  nvsl = FALSE;
596 
597  reduceStatsPrint(fullreduce, "nvsl", nvslnelims);
598 
599  if( SCIPgetTotalTime(scip) > timelimit )
600  break;
601  }
602 
603  if( bidecompparams && bidecompparams->depth < bidecompparams->maxdepth )
604  {
605  SCIPdebugMessage("go with depth %d \n", bidecompparams->depth);
606  SCIP_CALL( reduce_bidecomposition(scip, g, redbase, redbaseGetSolnode(redsollocal, redbase), wasDecomposed) );
607  SCIPdebugMessage("wasDecomposed=%d \n", *wasDecomposed);
608 
609  if( *wasDecomposed )
610  break;
611  }
612 
613  if( (inner_rounds > 0) && bred && nodereplacing )
614  {
615  SCIP_Real ub;
616  reduce_sollocalSetOffset(*redbaseGetOffsetPointer(redbase), redsollocal);
617  ub = reduce_sollocalGetUpperBound(redsollocal);
618  SCIP_CALL(reduce_bound(scip, g,
619  redbase->vnoi, redbase->nodearrreal, redbaseGetOffsetPointer(redbase), &ub,
620  redbase->heap, redbase->state, redbase->vbase, &brednelims));
621 
622  if( brednelims <= reductbound )
623  bred = FALSE;
624 
625  reduceStatsPrint(fullreduce, "bnd", brednelims);
626 
627  if( SCIPgetTotalTime(scip) > timelimit )
628  break;
629  }
630 
631  if( da )
632  {
633  const RPDA paramsda = {
634  .damode = (!redparameters->userec && inner_rounds > 0) ? STP_DAMODE_MEDIUM : STP_DAMODE_FAST,
635  .useRec = FALSE,
636  .useSlackPrune = FALSE, .extredMode = extred_none, .nodereplacing = nodereplacing };
637 
638  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal,
639  redbaseGetOffsetPointer(redbase), &danelims, randnumgen));
640 
641  if( danelims <= STP_RED_EXPENSIVEFACTOR * reductbound )
642  da = FALSE;
643 
644  reduceStatsPrint(fullreduce, "da", danelims);
645 
646  if( SCIPgetTotalTime(scip) > timelimit )
647  break;
648  }
649 
650  SCIP_CALL(reduce_unconnected(scip, g));
651  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &degtnelims, NULL));
652 
653  /* too few eliminations? */
654  if( (sdbiasnelims + danelims + sdnelims + bdknelims + nvslnelims + pathrepnelims +
655  brednelims + sdcnelims + sdstarnelims) <= STP_RED_GLBFACTOR * reductbound )
656  {
657  // at least one successful round and full reduce and no inner_restarts yet?
658  if( inner_rounds > 0 && fullreduce && inner_restarts == 0 )
659  {
660  inner_restarts++;
661  sd = TRUE;
662  sdstar = TRUE;
663  sdbiased = TRUE;
664  sdc = TRUE;
665  nvsl = nodereplacing;
666 
667  assert(extensive || sdcnelims == 0);
668 
669 #ifdef STP_PRINT_STATS
670  printf("RESTART inner reductions (restart number %d) \n", inner_restarts);
671 #endif
672  }
673  else
674  {
675  rerun = FALSE;
676  }
677  }
678 
679  if( extensive && (sdbiasnelims + danelims + sdnelims + bdknelims + nvslnelims + brednelims + sdcnelims + sdstarnelims) > 0 )
680  rerun = TRUE;
681 
682  inner_rounds++;
683  } /* inner reduction loop */
684 
685  return SCIP_OKAY;
686 }
687 
688 /** inner MWCS loop
689  * todo use REDBASE here instead of all the parameters and also allocate the memory locally! */
690 static
692  SCIP* scip, /**< SCIP data structure */
693  GRAPH* g, /**< graph data structure */
694  REDSOLLOCAL* redsollocal, /**< data structure for retaining primal solution */
695  PATH* vnoi, /**< Voronoi data structure */
696  SCIP_Real* nodearrreal, /**< nodes-sized array */
697  int* state, /**< shortest path array */
698  int* vbase, /**< voronoi base array */
699  int* nodearrint, /**< nodes-sized array */
700  STP_Bool* nodearrchar, /**< nodes-sized array */
701  SCIP_Real* fixed, /**< pointer to store the offset value */
702  STP_Bool dualascent, /**< do dual-ascent reduction? */
703  STP_Bool bred, /**< do bound-based reduction? */
704  int redbound, /**< minimal number of edges to be eliminated in order to reiterate reductions */
705  SCIP_Bool userec, /**< use recombination heuristic? */
706  SCIP_Real prizesum /**< prize sum */
707  )
708 {
709  int* solnode = reduce_sollocalGetSolnode(redsollocal);
710  SCIP_Real timelimit;
711  STP_Bool da = dualascent;
712  STP_Bool ans = TRUE;
713  STP_Bool nnp = TRUE;
714  STP_Bool npv = TRUE;
715  STP_Bool rerun = TRUE;
716  STP_Bool ansad = TRUE;
717  STP_Bool ansad2 = TRUE;
718  STP_Bool chain2 = TRUE;
719  STP_Bool extensive = STP_RED_EXTENSIVE;
720  SCIP_RANDNUMGEN* randnumgen;
721 
722  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
723  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
724 
725  /* inner loop */
726  for( int rounds = 0; rounds < STP_RED_MAXNINNROUNDS && !SCIPisStopped(scip) && rerun; rounds++ )
727  {
728  int daelims = 0;
729  int anselims = 0;
730  int nnpelims = 0;
731  int npvelims = 0;
732  int bredelims = 0;
733  int ansadelims = 0;
734  int ansad2elims = 0;
735  int chain2elims = 0;
736  int degelims = 0;
737 
738  if( SCIPgetTotalTime(scip) > timelimit )
739  break;
740 
741  if( ans || extensive )
742  {
743  SCIP_CALL( reduce_ans(scip, g, &anselims) );
744  if( anselims <= redbound )
745  ans = FALSE;
746  }
747 
748  if( ansad || extensive )
749  {
750  SCIP_CALL( reduce_ansAdv(scip, g, &ansadelims, FALSE) );
751  if( ansadelims <= redbound )
752  ansad = FALSE;
753 
754  SCIPdebugMessage("ans advanced deleted: %d \n", ansadelims);
755  }
756 
757  if( ans || ansad || nnp || npv || extensive )
758  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
759 
760  if( (da || (dualascent && extensive)) )
761  {
762  const RPDA paramsda = {
763  .damode = STP_DAMODE_FAST , .useRec = userec, .extredMode = extred_none,
764  .nodereplacing = FALSE, .useSlackPrune = FALSE,
765  .pcmw_solbasedda = (rounds > 0), .pcmw_useMultRoots = FALSE, .pcmw_markroots = FALSE, .pcmw_fastDa = (rounds == 0) };
766 
767  if( graph_pc_isRootedPcMw(g) )
768  {
769  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &daelims, randnumgen) );
770  }
771  else
772  {
773  reduce_sollocalSetOffset(*fixed, redsollocal);
774  SCIP_CALL( reduce_sollocalRebuildTry(scip, g, redsollocal) );
775 
776  SCIP_CALL( reduce_daPcMw(scip, g, &paramsda, redsollocal,
777  vnoi, nodearrreal, vbase, nodearrint, state, nodearrchar, &daelims, randnumgen, prizesum) );
778  }
779 
780  if( daelims <= 2 * redbound )
781  da = FALSE;
782  else
783  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
784 
785  SCIPdebugMessage("Dual-ascent eliminations: %d \n", daelims);
786  }
787 
788  if( nnp )
789  {
790  SCIP_CALL( reduce_nnp(scip, g, &nnpelims) );
791  if( nnpelims <= redbound )
792  nnp = FALSE;
793 
794  SCIPdebugMessage("nnp eliminations: %d \n", nnpelims);
795  }
796 
797  if( nnp || extensive )
798  {
799  SCIP_CALL(reduce_chain2(scip, g, vnoi, state, vbase, nodearrint, &chain2elims, 500));
800  if( chain2elims <= redbound )
801  chain2 = FALSE;
802 
803  SCIPdebugMessage("chain2 eliminations: %d \n", chain2elims);
804 
805  if( SCIPgetTotalTime(scip) > timelimit )
806  break;
807  }
808 
809  if( npv || extensive )
810  {
811  SCIP_CALL(reduce_npv(scip, g, vnoi, state, vbase, nodearrint, &npvelims, 400));
812  if( npvelims <= redbound )
813  npv = FALSE;
814 
815  SCIPdebugMessage("npv eliminations: %d \n", npvelims);
816  }
817 
818  if( chain2 || extensive )
819  {
820  SCIP_CALL(reduce_chain2(scip, g, vnoi, state, vbase, nodearrint, &chain2elims, 300));
821  if( chain2elims <= redbound )
822  chain2 = FALSE;
823 
824  SCIPdebugMessage("chain2 eliminations: %d \n", chain2elims);
825  }
826 
827  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
828 
829  if( ansad2 || extensive )
830  {
831  SCIP_CALL( reduce_ansAdv2(scip, g, &ansad2elims) );
832  if( ansad2elims <= redbound )
833  ansad2 = FALSE;
834  else
835  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &ansad2elims) );
836 
837  SCIPdebugMessage("ans advanced 2 eliminations: %d (da? %d ) \n", ansad2elims, da);
838  }
839 
840  if( bred )
841  {
842  SCIP_CALL( reduce_boundMw(scip, g, vnoi, fixed, nodearrint, state, vbase, NULL, &bredelims) );
843  if( bredelims <= redbound )
844  bred = FALSE;
845 
846  SCIPdebugMessage("reduce_bound eliminations: %d \n", bredelims);
847  }
848 
849  if( anselims + nnpelims + chain2elims + bredelims + npvelims + ansadelims + ansad2elims + daelims <= redbound )
850  rerun = FALSE;
851 
852 #ifdef SCIP_DEBUG
853  if( advanced && tryrmw )
854  {
855  int nreal, mreal;
856  graph_get_nVET(g, &nreal, &mreal, NULL);
857 
858  printf("round %d (of %d) nreal=%d mreal=%d\n", rounds, STP_RED_MAXNINNROUNDS, nreal, mreal);
859  }
860 #endif
861  } /* inner loop */
862 
863  SCIPfreeRandom(scip, &randnumgen);
864 
865  return SCIP_OKAY;
866 }
867 
868 
869 /** inner (R)PC loop
870  * todo use REDBASE here instead of all the parameters and also allocate the memory locally! */
871 static
873  SCIP* scip, /**< SCIP data structure */
874  GRAPH* g, /**< graph data structure */
875  REDSOLLOCAL* redsollocal, /**< solution store */
876  DHEAP* dheap, /**< heap data structure */
877  PATH* vnoi, /**< Voronoi data structure */
878  PATH* path, /**< path data structure */
879  SCIP_Real* nodearrreal, /**< nodes-sized array */
880  int* heap, /**< shortest path array */
881  int* state, /**< voronoi base array */
882  int* vbase, /**< nodes-sized array */
883  int* nodearrint, /**< node-sized array */
884  int* edgearrint, /**< edge-sized array */
885  int* nodearrint2, /**< nodes-sized array */
886  STP_Bool* nodearrchar, /**< nodes-sized array */
887  SCIP_Real* fixed, /**< offset */
888  SCIP_RANDNUMGEN* randnumgen, /**< random number generator */
889  SCIP_Real prizesum, /**< prize sum */
890  SCIP_Bool dualascent, /**< do dual-ascent reduction? */
891  SCIP_Bool bred, /**< do bound-based reduction? */
892  int reductbound, /**< minimal number of edges to be eliminated in order to reiterate reductions */
893  SCIP_Bool userec, /**< use recombination heuristic? */
894  SCIP_Bool nodereplacing, /**< should node replacement (by edges) be performed? */
895  SCIP_Bool usestrongreds, /**< allow strong reductions? */
896  SCIP_Bool beFast, /**< fast mode? */
897  int* ninnerelims /**< number of eliminations */
898  )
899 {
900  int* solnode = reduce_sollocalGetSolnode(redsollocal);
901  SCIP_Real timelimit;
902  const int reductbound_global = reductbound * STP_RED_GLBFACTOR;
903  SCIP_Bool dapaths = usestrongreds;
904  SCIP_Bool da = dualascent;
905  SCIP_Bool sd = usestrongreds;
906  SCIP_Bool sdw = usestrongreds;
907  SCIP_Bool sdstar = TRUE;
908  SCIP_Bool bd3 = nodereplacing;
909  SCIP_Bool nvsl = TRUE;
910  SCIP_Bool pathrep = usestrongreds;
911  SCIP_Bool rerun = TRUE;
912  const SCIP_Bool verbose = FALSE && dualascent && userec && nodereplacing;
913  const STP_Bool extensive = STP_RED_EXTENSIVE;
914 
915  *ninnerelims = 0;
916 
917  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
918 
919  /* inner loop! */
920  for( int rounds = 0; rounds < STP_RED_MAXNINNROUNDS && !SCIPisStopped(scip) && rerun; rounds++ )
921  {
922  int dapathelims = 0;
923  int danelims = 0;
924  int sdnelims = 0;
925  int sdcnelims = 0;
926  int bd3nelims = 0;
927  int nvslnelims = 0;
928  int sdwnelims = 0;
929  int sdstarnelims = 0;
930  int brednelims = 0;
931  int degnelims = 0;
932  int nelims = 0;
933  int pathrepnelims = 0;
934 
935 
936  if( SCIPgetTotalTime(scip) > timelimit )
937  break;
938 
939  if( sdstar || extensive )
940  {
941  int sdstarpcnelims = 0;
942 
943  SCIP_CALL( reduce_sdStarBiased(scip, getWorkLimitsPc(g, rounds, beFast, pc_sdstar), usestrongreds, g, &sdstarnelims));
944  if( verbose ) printf("sdstarnelims %d \n", sdstarnelims);
945 
946  if( sdstarnelims <= reductbound )
947  {
948  sdstar = FALSE;
949  }
950  else
951  {
952  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &nelims, &degnelims, solnode) );
953  SCIP_CALL( reduce_sdStarPc2(scip, getWorkLimitsPc(g, rounds, beFast, pc_sdstar2), usestrongreds, g, nodearrreal, nodearrint, nodearrint2, nodearrchar, dheap, &sdstarpcnelims));
954  if( verbose ) printf("sdstarpcnelims %d \n", sdstarpcnelims);
955 
956  sdstarnelims += sdstarpcnelims;
957 
958  }
959  }
960 
961  if( (graph_pc_isRootedPcMw(g) && (dapaths || extensive)) )
962  {
963  reduce_dapaths(scip, g, fixed, &dapathelims);
964 
965  if( dapathelims <= reductbound )
966  dapaths = FALSE;
967  }
968 
969  if( sd || extensive )
970  {
971  SCIP_CALL( execPc_SD(scip, g, vnoi, heap, state, vbase, nodearrint, nodearrint2, &sdnelims,
972  reductbound, verbose, &sd) );
973  }
974 
975  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &nelims, &degnelims, solnode) );
976 
977  if( sdw || extensive )
978  {
979  SCIP_CALL( reduce_sdWalkTriangle(scip, getWorkLimitsPc(g, rounds, beFast, pc_sdw1), usestrongreds, g, nodearrint, nodearrreal, vbase, nodearrchar, dheap, &sdwnelims));
980 
981  if( verbose ) printf("SDw: %d\n", sdwnelims);
982 
983  if( sdwnelims <= reductbound )
984  {
985  sdw = FALSE;
986  }
987  else
988  {
989  int sdwnelims2 = 0;
990  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &nelims, &degnelims, solnode) );
991  SCIP_CALL( reduce_sdWalkExt(scip, getWorkLimitsPc(g, rounds, beFast, pc_sdw2), usestrongreds, g, nodearrreal, heap, state, vbase, nodearrchar, &sdwnelims2) );
992  sdwnelims += sdwnelims2;
993 
994  if( verbose ) printf("SDwExt: %d\n", sdwnelims2);
995  }
996  }
997 
998  if( SCIPgetTotalTime(scip) > timelimit )
999  break;
1000 
1001  if( bd3 || extensive )
1002  {
1003  SCIP_CALL( execPc_BDk(scip, g, vnoi, path, heap, state, vbase, nodearrint, nodearrint2,
1004  &bd3nelims, getWorkLimitsPc(g, rounds, beFast, pc_bd3), fixed, reductbound, verbose, &bd3) );
1005  }
1006 
1007  if( nvsl || extensive )
1008  {
1009  SCIP_CALL( execPc_NVSL(scip, usestrongreds, g, vnoi, nodearrreal, fixed, edgearrint, vbase,
1010  nodearrint, nodearrint2, solnode, nodearrchar, &nvslnelims, reductbound, verbose, &nvsl) );
1011  }
1012 
1013  if( bred || extensive)
1014  {
1015  SCIP_CALL( execPc_BND(scip, g, vnoi, nodearrreal, fixed, heap, state, vbase, &brednelims, reductbound, verbose, &bred) );
1016  }
1017 
1018  if( SCIPgetTotalTime(scip) > timelimit )
1019  break;
1020 
1021  if( da || (dualascent && extensive) )
1022  {
1023  RPDA paramsda = { .damode = STP_DAMODE_FAST, .useRec = userec, .extredMode = extred_none, .nodereplacing = nodereplacing, .useSlackPrune = FALSE,
1024  .pcmw_solbasedda = !beFast,
1025  .pcmw_useMultRoots = FALSE, .pcmw_markroots = FALSE,
1026  .pcmw_fastDa = beFast };
1027 
1028  if( g->edges > STP_RED_EDGELIMIT_HUGE )
1029  {
1030  paramsda.pcmw_solbasedda = TRUE;
1031  paramsda.pcmw_fastDa = (beFast && rounds == 0);
1032  }
1033 
1034  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &nelims, &degnelims, solnode) );
1035 
1036  if( g->stp_type == STP_RPCSPG )
1037  {
1038  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &danelims, randnumgen) );
1039  }
1040  else
1041  {
1042  reduce_sollocalSetOffset(*fixed, redsollocal);
1043  SCIP_CALL( reduce_sollocalRebuildTry(scip, g, redsollocal) );
1044 
1045  SCIP_CALL( reduce_daPcMw(scip, g, &paramsda, redsollocal, vnoi, nodearrreal, vbase, heap,
1046  state, nodearrchar, &danelims, randnumgen, prizesum) );
1047  }
1048 
1049  if( danelims <= reductbound )
1050  da = FALSE;
1051 
1052  if( verbose ) printf("daX: %d \n", danelims);
1053  }
1054 
1055  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &nelims, &degnelims, solnode) );
1056 
1057  if( pathrep )
1058  {
1059  SCIP_CALL( reduce_pathreplace(scip, g, &pathrepnelims) );
1060 
1061  if( pathrepnelims <= reductbound )
1062  pathrep = FALSE;
1063  if( verbose ) printf("pathrep: %d\n", pathrepnelims);
1064  }
1065 
1066  nelims = degnelims + sdnelims + sdcnelims + bd3nelims + danelims + brednelims + nvslnelims
1067  + sdwnelims + sdstarnelims + dapathelims + pathrepnelims;
1068  *ninnerelims += nelims;
1069 
1070  if( nelims <= reductbound_global )
1071  rerun = FALSE;
1072  } /* inner reduction loop */
1073 
1074  return SCIP_OKAY;
1075 }
1076 
1077 
1078 
1079 
1080 /*
1081  * Interface methods
1082  */
1083 
1084 
1085 
1086 /** gets reduction bound */
1088  const GRAPH* g, /**< graph data structure */
1089  int lowerbound /**< lower bound on number of reductions (>= 1) */
1090  )
1091 {
1092  int min;
1093  const int nedges = graph_get_nEdges(g);
1094  const int nnodes = graph_get_nNodes(g);
1095 
1096  assert(lowerbound >= 1);
1097 
1098  if( graph_typeIsSpgLike(g) )
1099  {
1100  min = MAX(nedges / 1000, lowerbound);
1101  }
1102  else if( graph_pc_isPc(g) )
1103  {
1104  min = MAX(nnodes / 500, lowerbound);
1105  }
1106  else if( graph_pc_isMw(g) )
1107  {
1108  min = MAX(nnodes / 500, lowerbound);
1109  }
1110  else
1111  {
1112  min = MAX(nnodes / 1000, lowerbound);
1113  }
1114 
1115  return min;
1116 }
1117 
1118 
1119 /** initializes */
1121  SCIP* scip, /**< SCIP data structure */
1122  const GRAPH* g, /**< graph data structure */
1123  REDBASE** redbase /**< base */
1124  )
1125 {
1126  REDBASE* red;
1127  const int nnodes = graph_get_nNodes(g);
1128  const int nedges = graph_get_nEdges(g);
1129 
1130  assert(scip);
1131 
1132  SCIP_CALL( SCIPallocMemory(scip, redbase) );
1133  red = *redbase;
1134 
1135  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->edgearrint), nedges) );
1136  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->nodearrchar), nnodes) );
1137  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->heap), nnodes + 1) );
1138  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->state), 4 * nnodes) );
1139  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->nodearrreal), nnodes) );
1140  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->vbase), 4 * nnodes) );
1141  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->nodearrint), nnodes) );
1142  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->nodearrint2), nnodes) );
1143  SCIP_CALL( SCIPallocMemoryArray(scip, &(red->vnoi), 4 * nnodes) );
1144  red->path = NULL;
1145 
1146  red->redparameters = NULL;
1147  red->solnode = NULL;
1148  red->redsol = NULL;
1149 
1150  return SCIP_OKAY;
1151 }
1152 
1153 
1154 /** frees */
1156  SCIP* scip, /**< SCIP data structure */
1157  REDBASE** redbase /**< base */
1158  )
1159 {
1160  REDBASE* red;
1161  assert(scip && redbase);
1162  assert(*redbase);
1163 
1164  red = *redbase;
1165 
1166  SCIPfreeMemoryArray(scip, &(red->vnoi));
1167  SCIPfreeMemoryArray(scip, &(red->nodearrint2));
1168  SCIPfreeMemoryArray(scip, &(red->nodearrint));
1169  SCIPfreeMemoryArray(scip, &(red->vbase));
1170  SCIPfreeMemoryArray(scip, &(red->nodearrreal));
1171  SCIPfreeMemoryArray(scip, &(red->state));
1172  SCIPfreeMemoryArray(scip, &(red->heap));
1173  SCIPfreeMemoryArray(scip, &(red->nodearrchar));
1174  SCIPfreeMemoryArray(scip, &(red->edgearrint));
1175 
1176  SCIPfreeMemory(scip, redbase);
1177 }
1178 
1179 
1180 /** basic reduction package for the STP */
1182  SCIP* scip, /**< SCIP data structure */
1183  GRAPH* g, /**< graph data structure */
1184  REDSOL* redsol, /**< primal solution container */
1185  int minelims, /**< minimal number of edges to be eliminated in order to reiterate reductions */
1186  SCIP_Bool dualascent, /**< perform dual-ascent reductions? */
1187  SCIP_Bool nodereplacing, /**< should node replacement (by edges) be performed? */
1188  SCIP_Bool userec, /**< use recombination heuristic? */
1189  SCIP_Bool usestrongreds /**< allow strong reductions?
1190  NOTE: needed for propagation, because arcs might have been fixed to 0 */
1191  )
1192 {
1193  PATH* vnoi;
1194  PATH* path;
1195  SCIP_Real* nodearrreal;
1196  int* heap;
1197  int* state;
1198  int* vbase;
1199  int* nodearrint;
1200  int* edgearrint;
1201  int* nodearrint2;
1202  STP_Bool* nodearrchar;
1203  const int nnodes = graph_get_nNodes(g);
1204  const int nedges = graph_get_nEdges(g);
1205  const int nterms = graph_get_nTerms(g);
1206  const int reductbound = reduce_getMinNreductions(g, minelims);
1207  SCIP_Bool bred = FALSE;
1208 
1209  assert(scip && redsol);
1210 
1211  if( SCIPisLE(scip, (SCIP_Real) nterms / (SCIP_Real) nnodes, STP_BND_THRESHOLD) )
1212  bred = TRUE;
1213 
1214  /* allocate memory */
1215  SCIP_CALL( SCIPallocBufferArray(scip, &edgearrint, nedges) );
1216  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrchar, nnodes) );
1217  SCIP_CALL( SCIPallocBufferArray(scip, &heap, nnodes + 1) );
1218  SCIP_CALL( SCIPallocBufferArray(scip, &state, 4 * nnodes) );
1219  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrreal, nnodes) );
1220  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, 4 * nnodes) );
1221  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint, nnodes) );
1222  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint2, nnodes) );
1223  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, 4 * nnodes) );
1224  path = NULL;
1225 
1226  /* reduction loop */
1227  {
1228  RPARAMS parameters = { .dualascent = dualascent, .boundreduce = bred, .nodereplacing = nodereplacing,
1229  .reductbound_min = minelims,
1230  .reductbound = reductbound, .userec = userec, .fullreduce = (dualascent && userec),
1231  .usestrongreds = usestrongreds };
1232  BIDECPARAMS decparameters = { .depth = 0, .maxdepth = 3, .newLevelStarted = FALSE };
1233  REDBASE redbase = { .redparameters = &parameters, .bidecompparams = &decparameters,
1234  .solnode = NULL, .redsol = redsol,
1235  .vnoi = vnoi, .path = path, .heap = heap,
1236  .nodearrreal = nodearrreal,
1237  .state = state, .vbase = vbase, .nodearrint = nodearrint,
1238  .edgearrint = edgearrint, .nodearrint2 = nodearrint2, .nodearrchar = nodearrchar };
1239 
1240  SCIP_CALL( reduce_redLoopStp(scip, g, &redbase) );
1241  }
1242 
1243  SCIPdebugMessage("Reduction Level 1: Fixed Cost = %.12e\n", reduce_solGetOffset(redsol));
1244 
1245  /* free memory */
1246  SCIPfreeBufferArray(scip, &vnoi);
1247  SCIPfreeBufferArray(scip, &nodearrint2);
1248  SCIPfreeBufferArray(scip, &nodearrint);
1249  SCIPfreeBufferArray(scip, &vbase);
1250  SCIPfreeBufferArray(scip, &nodearrreal);
1251  SCIPfreeBufferArray(scip, &state);
1252  SCIPfreeBufferArray(scip, &heap);
1253  SCIPfreeBufferArray(scip, &nodearrchar);
1254  SCIPfreeBufferArray(scip, &edgearrint);
1255 
1256  return SCIP_OKAY;
1257 }
1258 
1259 /** basic reduction package for the (R)PCSTP */
1261  SCIP* scip, /**< SCIP data structure */
1262  REDSOL* redsol, /**< solution storage */
1263  GRAPH* g, /**< graph data structure */
1264  int minelims, /**< minimal number of edges to be eliminated in order to reiterate reductions */
1265  SCIP_Bool advanced, /**< perform advanced (e.g. dual ascent) reductions? */
1266  SCIP_Bool userec, /**< use recombination heuristic? */
1267  SCIP_Bool nodereplacing, /**< should node replacement (by edges) be performed? */
1268  SCIP_Bool usestrongreds /**< allow strong reductions?
1269  NOTE: needed for propagation, because arcs might have been fixed to 0 */
1270  )
1271 {
1272  PATH* vnoi;
1273  PATH* path;
1274  SCIP_Real* nodearrreal;
1275  SCIP_Real timelimit;
1276  int* heap;
1277  int* state;
1278  int* vbase;
1279  int* nodearrint;
1280  int* edgearrint;
1281  int* nodearrint2;
1282  int nnodes;
1283  int nterms;
1284  int nedges;
1285  int reductbound;
1286  STP_Bool* nodearrchar;
1287  SCIP_Bool bred = FALSE;
1288 
1289  assert(scip != NULL);
1290  assert(g != NULL);
1291  assert(minelims >= 0);
1292 
1293  nterms = g->terms;
1294  nnodes = g->knots;
1295  nedges = g->edges;
1296 
1297  /* get timelimit parameter*/
1298  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1299 
1300  /* allocate memory */
1301  SCIP_CALL( SCIPallocBufferArray(scip, &heap, nnodes + 1) );
1302  SCIP_CALL( SCIPallocBufferArray(scip, &state, 4 * nnodes) );
1303  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrreal, nnodes + 2) );
1304  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, 4 * nnodes) );
1305  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, 4 * nnodes) );
1306  SCIP_CALL( SCIPallocBufferArray(scip, &path, nnodes + 1) );
1307  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint, nnodes + 1) );
1308  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint2, nnodes + 1) );
1309  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrchar, nnodes + 1) );
1310  SCIP_CALL( SCIPallocBufferArray(scip, &edgearrint, nedges) );
1311 
1312  if( SCIPisLE(scip, (SCIP_Real) nterms / (SCIP_Real) nnodes, STP_BND_THRESHOLD) )
1313  {
1314  bred = TRUE;
1315 
1316  if( SCIPisGT(scip, (SCIP_Real) (nterms + graph_pc_nNonLeafTerms(g)) / (SCIP_Real) nnodes, STP_BND_THRESHOLD ) )
1317  bred = FALSE;
1318  }
1319 
1320  /* define minimal number of edge/node eliminations for a reduction test to be continued */
1321  reductbound = reduce_getMinNreductions(g, minelims);
1322 
1323  /* reduction loop */
1324  SCIP_CALL( reduce_redLoopPc(scip, redsol, g, vnoi, path, nodearrreal, heap, state,
1325  vbase, nodearrint, edgearrint, nodearrint2, nodearrchar,
1326  advanced, bred, userec && advanced, reductbound, userec, nodereplacing, usestrongreds) );
1327 
1328  /* free memory */
1329  SCIPfreeBufferArray(scip, &edgearrint);
1330  SCIPfreeBufferArray(scip, &nodearrchar);
1331  SCIPfreeBufferArray(scip, &nodearrint2);
1332  SCIPfreeBufferArray(scip, &nodearrint);
1333  SCIPfreeBufferArray(scip, &path);
1334  SCIPfreeBufferArray(scip, &vnoi);
1335  SCIPfreeBufferArray(scip, &vbase);
1336  SCIPfreeBufferArray(scip, &nodearrreal);
1337  SCIPfreeBufferArray(scip, &state);
1338  SCIPfreeBufferArray(scip, &heap);
1339 
1340  return SCIP_OKAY;
1341 }
1342 
1343 /** reduction package for the MWCSP */
1345  SCIP* scip, /**< SCIP data structure */
1346  REDSOL* redsol, /**< solution storage */
1347  GRAPH* g, /**< graph data structure */
1348  int minelims, /**< minimal number of edges to be eliminated in order to reiterate reductions */
1349  SCIP_Bool advanced, /**< perform advanced reductions? */
1350  SCIP_Bool userec, /**< use recombination heuristic? */
1351  SCIP_Bool usestrongreds /**< allow strong reductions?
1352  NOTE: needed for propagation, because arcs might have been fixed to 0 */
1353  )
1354 {
1355  PATH* vnoi;
1356  SCIP_Real* nodearrreal;
1357  int* state;
1358  int* vbase;
1359  int* nodearrint;
1360  const int nnodes = graph_get_nNodes(g);
1361  const int nterms = graph_get_nTerms(g);
1362  const int redbound = reduce_getMinNreductions(g, minelims);
1363  STP_Bool* nodearrchar;
1364  STP_Bool bred = FALSE;
1365 
1366  assert(scip != NULL);
1367  assert(redsol != NULL);
1368 
1369  if( SCIPisLE(scip, (SCIP_Real) nterms / (SCIP_Real) nnodes, 0.1) )
1370  bred = TRUE;
1371 
1372  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint, nnodes + 1) );
1373  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrchar, nnodes + 1) );
1374  SCIP_CALL( SCIPallocBufferArray(scip, &state, 4 * nnodes) );
1375  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, 4 * nnodes) );
1376  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, 4 * nnodes) );
1377 
1378  if( bred || advanced )
1379  {
1380  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrreal, nnodes + 2) );
1381  }
1382  else
1383  {
1384  nodearrreal = NULL;
1385  }
1386 
1387  /* reduction loop */
1388  SCIP_CALL( reduce_redLoopMw(scip, redsol, g, vnoi, nodearrreal, state,
1389  vbase, nodearrint, nodearrchar, advanced, bred, advanced, redbound, userec, usestrongreds) );
1390 
1391  /* free memory */
1392  SCIPfreeBufferArrayNull(scip, &nodearrreal);
1393  SCIPfreeBufferArray(scip, &vnoi);
1394  SCIPfreeBufferArray(scip, &vbase);
1395  SCIPfreeBufferArray(scip, &state);
1396  SCIPfreeBufferArray(scip, &nodearrchar);
1397  SCIPfreeBufferArray(scip, &nodearrint);
1398 
1399  return SCIP_OKAY;
1400 }
1401 
1402 /** basic reduction package for the HCDSTP */
1404  SCIP* scip, /**< SCIP data structure */
1405  GRAPH* g, /**< graph data structure */
1406  SCIP_Real* fixed, /**< pointer to store the offset value */
1407  int minelims /**< minimal number of edges to be eliminated in order to reiterate reductions */
1408  )
1409 {
1410  REDSOLLOCAL* redsollocal;
1411  SCIP_RANDNUMGEN* randnumgen;
1412  PATH* vnoi;
1413  SCIP_Real* cost;
1414  SCIP_Real* radius;
1415  SCIP_Real* costrev;
1416  SCIP_Real timelimit;
1417  int* heap;
1418  int* state;
1419  int* vbase;
1420  int* pathedge;
1421  STP_Bool* nodearrchar;
1422  int redbound;
1423  int degnelims;
1424  STP_Bool bred = FALSE; // todo reduction method is not correct! ignores direction
1425  STP_Bool hbred = TRUE;
1426  STP_Bool rbred = TRUE;
1427  STP_Bool rcbred = TRUE;
1428  STP_Bool da = TRUE;
1429  STP_Bool dahop = TRUE;
1430  int nrounds = 0;
1431  const int nnodes = graph_get_nNodes(g);
1432  const int nedges = graph_get_nEdges(g);
1433 
1434  assert(scip != NULL);
1435  assert(g != NULL);
1436  assert(minelims >= 0);
1437 
1438  degnelims = 0;
1439  redbound = MAX(g->knots / 1000, minelims);
1440  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1441 
1442  /* allocate memory */
1443  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrchar, nnodes) );
1444  SCIP_CALL( SCIPallocBufferArray(scip, &heap, nnodes + 1) );
1445  SCIP_CALL( SCIPallocBufferArray(scip, &state, 3 * nnodes) );
1446  SCIP_CALL( SCIPallocBufferArray(scip, &cost, nedges) );
1447  SCIP_CALL( SCIPallocBufferArray(scip, &radius, nnodes) );
1448  SCIP_CALL( SCIPallocBufferArray(scip, &costrev, nedges) );
1449  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, 3 * nnodes) );
1450  SCIP_CALL( SCIPallocBufferArray(scip, &pathedge, nnodes) );
1451  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, 3 * nnodes) );
1452 
1453  SCIP_CALL( reduce_simple_hc(scip, g, fixed, &degnelims) );
1455 
1456  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1457  SCIP_CALL( reduce_sollocalInit(scip, g, &redsollocal) );
1458 
1459  while( (da || bred || hbred || rbred || rcbred || dahop ) && !SCIPisStopped(scip) )
1460  {
1461  int ndahopelims = 0;
1462  int danelims = 0;
1463  int brednelims = 0;
1464  int hbrednelims = 0;
1465  int hcrnelims = 0;
1466  int hcrcnelims = 0;
1467 
1468  if( SCIPgetTotalTime(scip) > timelimit )
1469  break;
1470 
1471 
1472  if( rbred )
1473  {
1474  SCIP_CALL( reduce_boundHopR(scip, g, vnoi, cost, costrev, radius, heap, state, vbase, &hcrnelims, pathedge) );
1475  if( hcrnelims <= redbound )
1476  rbred = FALSE;
1477  }
1478 
1479  if( rcbred )
1480  {
1481  SCIP_CALL( reduce_boundHopRc(scip, g, vnoi, cost, costrev, radius, -1.0, heap, state, vbase, &hcrcnelims, pathedge, FALSE) );
1482  if( hcrcnelims <= redbound )
1483  rcbred = FALSE;
1484  }
1485 
1487 
1488  if( bred )
1489  {
1490  SCIP_Real upperbound = -1.0;
1491 
1492  SCIP_CALL( reduce_bound(scip, g, vnoi, radius, fixed, &upperbound, heap, state, vbase, &brednelims) );
1493  if( brednelims <= redbound )
1494  bred = FALSE;
1495  }
1496 
1497  if( SCIPgetTotalTime(scip) > timelimit )
1498  break;
1499 
1500  if( hbred )
1501  {
1502  SCIP_CALL( reduce_boundHop(scip, g, vnoi, cost, radius, costrev, heap, state, vbase, &hbrednelims) );
1503  if( hbrednelims <= redbound )
1504  hbred = FALSE;
1505  if( SCIPgetTotalTime(scip) > timelimit )
1506  break;
1507  }
1508 
1509  if( da )
1510  {
1511  const RPDA paramsda = { .damode = STP_DAMODE_FAST , .useSlackPrune = FALSE, .useRec = FALSE, .extredMode = extred_none, .nodereplacing = FALSE};
1512 
1513  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &danelims, randnumgen) );
1514  if( danelims <= redbound )
1515  da = FALSE;
1516  }
1517 
1518  if( SCIPgetTotalTime(scip) > timelimit )
1519  break;
1520 
1521  if( dahop )
1522  {
1523  SCIP_CALL( reduce_boundHopDa(scip, g, &ndahopelims, randnumgen) );
1524 
1525  if( ndahopelims <= redbound )
1526  dahop = FALSE;
1527  }
1528 
1529 
1530  if( da && nrounds == 0 )
1531  {
1532  rbred = TRUE;
1533  hbred = TRUE;
1534  rcbred = TRUE;
1535  }
1536 
1537  nrounds++;
1538  }
1539 
1540  reduce_sollocalFree(scip, &redsollocal);
1541  SCIPfreeRandom(scip, &randnumgen);
1542 
1543  SCIPfreeBufferArray(scip, &vnoi);
1544  SCIPfreeBufferArray(scip, &pathedge);
1545  SCIPfreeBufferArray(scip, &vbase);
1546  SCIPfreeBufferArray(scip, &costrev);
1547  SCIPfreeBufferArray(scip, &radius);
1548  SCIPfreeBufferArray(scip, &cost);
1549  SCIPfreeBufferArray(scip, &state);
1550  SCIPfreeBufferArray(scip, &heap);
1551  SCIPfreeBufferArray(scip, &nodearrchar);
1552 
1553  return SCIP_OKAY;
1554 }
1555 
1556 /** basic reduction package for the SAP */
1558  SCIP* scip, /**< SCIP data structure */
1559  GRAPH* g, /**< graph data structure */
1560  SCIP_Bool dualascent, /**< perform dual-ascent reductions? */
1561  SCIP_Real* fixed, /**< pointer to store the offset value */
1562  int minelims /**< minimal number of edges to be eliminated in order to reiterate reductions */
1563  )
1564 {
1565  PATH* vnoi;
1566  PATH* path;
1567  SCIP_Real timelimit;
1568  SCIP_Real* nodearrreal;
1569  int* heap;
1570  int* state;
1571  int* vbase;
1572  int* nodearrint;
1573  int* nodearrint2;
1574  int degtnelims;
1575  int redbound;
1576  const int nnodes = graph_get_nNodes(g);
1577  const int nedges = graph_get_nEdges(g);
1578  STP_Bool da = dualascent; // && (g->stp_type != STP_NWPTSPG);
1579  STP_Bool sd = TRUE;
1580  STP_Bool rpt = TRUE;
1581  SCIP_RANDNUMGEN* randnumgen;
1582 
1583  assert(!graph_typeIsUndirected(g));
1584 
1585  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1586 
1587  redbound = MAX(nnodes / 1000, minelims);
1588  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1589 
1590  /* allocate memory */
1591  SCIP_CALL( SCIPallocBufferArray(scip, &heap, nnodes + 1) );
1592  SCIP_CALL( SCIPallocBufferArray(scip, &state, nnodes) );
1593  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrreal, nnodes) );
1594  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, nnodes) );
1595  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint, nnodes) );
1596  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint2, nnodes) );
1597  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, nnodes) );
1598  SCIP_CALL( SCIPallocBufferArray(scip, &path, nnodes) );
1599 
1600  if( g->stp_type == STP_SAP )
1601  {
1602  /* todo change .stp file format for SAP! */
1603  for( int e = 0; e < nedges; e++ )
1604  {
1605  if( EQ(g->cost[e], 20000.0) )
1606  g->cost[e] = FARAWAY;
1607  }
1608  }
1609 
1610  SCIP_CALL( reduce_simple_sap(scip, g, fixed, &degtnelims) );
1611 
1612  /* main loop */
1613  while( (sd || rpt || da) && !SCIPisStopped(scip) )
1614  {
1615  int danelims = 0;
1616  int sdnelims = 0;
1617  int rptnelims = 0;
1618 
1619  if( SCIPgetTotalTime(scip) > timelimit )
1620  break;
1621 
1622  if( sd )
1623  {
1624  SCIP_CALL( reduce_sdspSap(scip, g, vnoi, path, heap, state, vbase, nodearrint, nodearrint2, &sdnelims, 300) );
1625 
1626  if( sdnelims <= redbound )
1627  sd = FALSE;
1628  }
1629 
1630  if( rpt )
1631  {
1632  SCIP_CALL( reduce_rpt(scip, g, fixed, &rptnelims) );
1633 
1634  if( rptnelims <= redbound )
1635  rpt = FALSE;
1636  }
1637 
1638  SCIP_CALL( reduce_simple_sap(scip, g, fixed, &degtnelims) );
1639 
1640  if( da )
1641  {
1642  const RPDA paramsda = { .damode = STP_DAMODE_FAST, .useSlackPrune = FALSE, .useRec = FALSE, .extredMode = extred_none, .nodereplacing = FALSE};
1643  SCIP_CALL( reduce_da(scip, g, &paramsda, NULL, fixed, &danelims, randnumgen) );
1644 
1645  if( danelims <= 2 * redbound )
1646  da = FALSE;
1647  }
1648  }
1649 
1650  SCIP_CALL( reduce_simple_sap(scip, g, fixed, &degtnelims) );
1651 
1652  SCIPfreeBufferArray(scip, &path);
1653  SCIPfreeBufferArray(scip, &vnoi);
1654  SCIPfreeBufferArray(scip, &nodearrint2);
1655  SCIPfreeBufferArray(scip, &nodearrint);
1656  SCIPfreeBufferArray(scip, &vbase);
1657  SCIPfreeBufferArray(scip, &nodearrreal);
1658  SCIPfreeBufferArray(scip, &state);
1659  SCIPfreeBufferArray(scip, &heap);
1660 
1661  /* free random number generator */
1662  SCIPfreeRandom(scip, &randnumgen);
1663 
1664  return SCIP_OKAY;
1665 }
1666 
1667 
1668 /** reduce node-weighted Steiner tree problem */
1670  SCIP* scip, /**< SCIP data structure */
1671  GRAPH* g, /**< graph data structure */
1672  SCIP_Real* fixed, /**< pointer to store the offset value */
1673  int minelims /**< minimal number of edges to be eliminated in order to reiterate reductions */
1674  )
1675 {
1676  PATH *vnoi;
1677  SCIP_Real *nodearrreal;
1678  SCIP_Real timelimit;
1679  int *vbase;
1680  int *nodearrint2;
1681  int nnodes;
1682  int redbound;
1683 
1684  STP_Bool da = TRUE;
1685  SCIP_RANDNUMGEN* randnumgen;
1686 
1687  /* create random number generator */
1688  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1689 
1690  nnodes = g->knots;
1691 
1692  redbound = MAX(nnodes / 1000, minelims);
1693  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1694 
1695  /* allocate memory */
1696  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrreal, nnodes) );
1697  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, nnodes) );
1698  SCIP_CALL( SCIPallocBufferArray(scip, &nodearrint2, nnodes) );
1699  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, nnodes) );
1700 
1701  while( (da) && !SCIPisStopped(scip) )
1702  {
1703  int danelims = 0;
1704  const RPDA paramsda = { .damode = STP_DAMODE_FAST , .useSlackPrune = FALSE, .useRec = FALSE, .extredMode = extred_none, .nodereplacing = FALSE};
1705 
1706  if( SCIPgetTotalTime(scip) > timelimit )
1707  break;
1708 
1709  SCIP_CALL( reduce_da(scip, g, &paramsda, NULL, fixed, &danelims, randnumgen) );
1710 
1711  if( danelims <= 2 * redbound )
1712  da = FALSE;
1713  }
1714 
1715  SCIPfreeBufferArray(scip, &vnoi);
1716  SCIPfreeBufferArray(scip, &nodearrint2);
1717  SCIPfreeBufferArray(scip, &vbase);
1718  SCIPfreeBufferArray(scip, &nodearrreal);
1719 
1720  /* free random number generator */
1721  SCIPfreeRandom(scip, &randnumgen);
1722 
1723  return SCIP_OKAY;
1724 }
1725 
1726 
1727 
1728 /** reduce degree constrained Steiner tree problem */
1730  SCIP* scip, /**< SCIP data structure */
1731  GRAPH* g, /**< graph data structure */
1732  SCIP_Real* fixed, /**< pointer to store the offset value */
1733  int minelims /**< minimal number of edges to be eliminated in order to reiterate reductions */
1734  )
1735 {
1736  REDSOLLOCAL* redsollocal;
1737  SCIP_RANDNUMGEN* randnumgen;
1738  SCIP_Real timelimit;
1739  const int nnodes = graph_get_nNodes(g);
1740  const int redbound = MAX(nnodes / 1000, minelims);
1741  STP_Bool da = TRUE;
1742 
1743  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1744  SCIP_CALL( reduce_sollocalInit(scip, g, &redsollocal) );
1745 
1746  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1747 
1748  reduce_simple_dc(scip, g);
1749 
1750  while( (da) && !SCIPisStopped(scip) )
1751  {
1752  int danelims = 0;
1753  const RPDA paramsda = { .damode = STP_DAMODE_FAST , .useSlackPrune = FALSE, .useRec = FALSE, .extredMode = extred_none, .nodereplacing = FALSE};
1754 
1755  if( SCIPgetTotalTime(scip) > timelimit )
1756  break;
1757 
1758  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &danelims, randnumgen) );
1759 
1760  if( danelims <= 2 * redbound )
1761  da = FALSE;
1762  }
1763 
1764  reduce_sollocalFree(scip, &redsollocal);
1765  SCIPfreeRandom(scip, &randnumgen);
1766 
1767  return SCIP_OKAY;
1768 }
1769 
1770 /** MWCS loop */
1772  SCIP* scip, /**< SCIP data structure */
1773  REDSOL* redsol, /**< solution contained */
1774  GRAPH* g, /**< graph data structure */
1775  PATH* vnoi, /**< Voronoi data structure */
1776  SCIP_Real* nodearrreal, /**< nodes-sized array */
1777  int* state, /**< shortest path array */
1778  int* vbase, /**< voronoi base array */
1779  int* nodearrint, /**< nodes-sized array */
1780  STP_Bool* nodearrchar, /**< nodes-sized array */
1781  STP_Bool advanced, /**< do advanced reduction? */
1782  STP_Bool bred, /**< do bound-based reduction? */
1783  STP_Bool tryrmw, /**< try to convert problem to RMWCSP? Only possible if advanced = TRUE and userec = TRUE */
1784  int redbound, /**< minimal number of edges to be eliminated in order to reiterate reductions */
1785  SCIP_Bool userec, /**< use recombination heuristic? */
1786  SCIP_Bool usestrongreds /**< allow strong reductions? */
1787 
1788  )
1789 {
1790  SCIP_Real* fixed;
1791  int* solnode;
1792  REDSOLLOCAL* redsollocal;
1793  SCIP_Real da = advanced;
1794  SCIP_Real timelimit;
1795  int degelims;
1796  SCIP_Bool fullrerun = TRUE;
1797  STP_Bool extensive = STP_RED_EXTENSIVE;
1798  SCIP_RANDNUMGEN* randnumgen;
1799  SCIP_Real prizesum;
1800 
1801  assert(scip && g);
1802  assert(advanced || !tryrmw);
1803 
1804  tryrmw = tryrmw && userec;
1805 
1806  SCIP_CALL( reduce_solInitLocal(scip, g, redsol, &redsollocal));
1807  fixed = reduce_solGetOffsetPointer(redsol);
1808  solnode = reduce_sollocalGetSolnode(redsollocal);
1809 
1810  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1811  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1812  graph_pc_2org(scip, g);
1813  degelims = 0;
1814  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
1815  assert(graph_pc_term2edgeIsConsistent(scip, g));
1816 
1817  prizesum = graph_pc_getPosPrizeSum(scip, g);
1818 
1819  /* todo currently run a most two times */
1820  for( int rounds = 0; rounds < 3 && !SCIPisStopped(scip) && fullrerun; rounds++ )
1821  {
1822  const RPDA paramsda = { .damode = STP_DAMODE_FAST , .useRec = userec, .extredMode = extred_none, .nodereplacing = FALSE, .useSlackPrune = TRUE,
1823  .pcmw_solbasedda = TRUE, .pcmw_useMultRoots = (g->terms > STP_RED_MWTERMBOUND),
1824  .pcmw_markroots = tryrmw, .pcmw_fastDa = FALSE };
1825  int daelims = 0;
1826 
1827  fullrerun = FALSE;
1828 
1829  SCIP_CALL( redLoopInnerMw(scip, g, redsollocal, vnoi, nodearrreal, state,
1830  vbase, nodearrint, nodearrchar, fixed, da, bred, redbound, userec, prizesum) );
1831 
1832  if( advanced && g->terms > 2 )
1833  {
1834  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
1835 
1836  if( graph_pc_isRootedPcMw(g) )
1837  {
1838  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &daelims, randnumgen) );
1839  }
1840  else
1841  {
1842  reduce_sollocalSetOffset(*fixed, redsollocal);
1843  SCIP_CALL( reduce_sollocalRebuildTry(scip, g, redsollocal) );
1844 
1845  SCIP_CALL( reduce_daPcMw(scip, g, &paramsda, redsollocal, vnoi, nodearrreal, vbase,
1846  nodearrint, state, nodearrchar, &daelims, randnumgen, prizesum) );
1847  }
1848 
1849  da = FALSE;
1850  userec = FALSE;
1851  advanced = FALSE;
1852  fullrerun = (daelims >= redbound || (extensive && (daelims > 0)));
1853 
1854  if( !graph_pc_isRootedPcMw(g) && tryrmw && g->terms > 2 )
1855  {
1856 #ifdef WITH_UG
1857  SCIP_CALL( graph_transPcmw2rooted(scip, g, prizesum, FALSE) );
1858 #else
1859  SCIP_CALL( graph_transPcmw2rooted(scip, g, prizesum, TRUE) );
1860 #endif
1861  if( graph_pc_isRootedPcMw(g) )
1862  {
1863  if( fullrerun )
1864  {
1865  advanced = TRUE;
1866  da = TRUE;
1867  }
1868  }
1869  }
1870 
1871  if( fullrerun )
1872  {
1873  SCIPdebugMessage("Restarting reduction loop after %d rounds! (%d DA eliminations) \n\n ", rounds, daelims);
1874 
1875  if( extensive )
1876  advanced = TRUE;
1877  }
1878  }
1879 
1880  SCIP_CALL( reduce_simple_mw(scip, g, solnode, fixed, &degelims) );
1881  }
1882 
1883  /* go back to the extended graph */
1884  graph_pc_2trans(scip, g);
1885 
1886  SCIP_CALL( reduce_unconnected(scip, g) );
1887 
1888  SCIPfreeRandom(scip, &randnumgen);
1889 
1890  reduce_solFinalizeLocal(scip, g, redsol);
1891 
1892  return SCIP_OKAY;
1893 }
1894 
1895 /** (R)PC loop */
1897  SCIP* scip, /**< SCIP data structure */
1898  REDSOL* redsol, /**< solution contained */
1899  GRAPH* g, /**< graph data structure */
1900  PATH* vnoi, /**< Voronoi data structure */
1901  PATH* path, /**< path data structure */
1902  SCIP_Real* nodearrreal, /**< nodes-sized array */
1903  int* heap, /**< shortest path array */
1904  int* state, /**< voronoi base array */
1905  int* vbase, /**< nodes-sized array */
1906  int* nodearrint, /**< node-sized array */
1907  int* edgearrint, /**< edge-sized array */
1908  int* nodearrint2, /**< nodes-sized array */
1909  STP_Bool* nodearrchar, /**< nodes-sized array */
1910  SCIP_Bool dualascent, /**< do dual-ascent reduction? */
1911  SCIP_Bool bred, /**< do bound-based reduction? */
1912  SCIP_Bool tryrpc, /**< try to transform to rpc? */
1913  int reductbound, /**< minimal number of edges to be eliminated in order to reiterate reductions */
1914  SCIP_Bool userec, /**< use recombination heuristic? */
1915  SCIP_Bool nodereplacing, /**< should node replacement (by edges) be performed? */
1916  SCIP_Bool usestrongreds /**< allow strong reductions? */
1917  )
1918 {
1919  SCIP_Real* fixed;
1920  int* solnode;
1921  REDSOLLOCAL* redsollocal;
1922  DHEAP* dheap;
1923  SCIP_RANDNUMGEN* randnumgen;
1924  SCIP_Real timelimit;
1925  SCIP_Bool rerun = TRUE;
1926  SCIP_Bool advancedrun = dualascent;
1927  SCIP_Real prizesum;
1928  // todo remove
1929  const SCIP_Bool verbose = FALSE && dualascent && userec && nodereplacing;
1930  int degnelims;
1931  int nadvruns = 0;
1932  SCIP_Bool isSpg = FALSE;
1933  SCIP_Bool rpc = (g->stp_type == STP_RPCSPG);
1934  const SCIP_Bool trySpgTrans = userec; // todo extra parameter
1935  const int reductbound_global = reductbound * STP_RED_GLBFACTOR;
1936 
1937  SCIP_CALL( reduce_solInitLocal(scip, g, redsol, &redsollocal));
1938  fixed = reduce_solGetOffsetPointer(redsol);
1939  solnode = reduce_sollocalGetSolnode(redsollocal);
1940 
1941  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
1942  SCIP_CALL( graph_heap_create(scip, g->knots, NULL, NULL, &dheap) );
1943  graph_pc_2org(scip, g);
1944 
1945  assert(graph_pc_term2edgeIsConsistent(scip, g));
1946 
1947  SCIP_CALL( graph_pc_presolInit(scip, g) );
1948 
1949  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &degnelims, NULL, solnode) );
1950  if( verbose ) printf("initial degnelims: %d \n", degnelims);
1951 
1952  prizesum = graph_pc_getPosPrizeSum(scip, g);
1953  assert(prizesum < FARAWAY);
1954 
1955  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1956 
1957  /* main reduction loop */
1958  for( int outterrounds = 0; outterrounds < STP_RED_MAXNOUTROUNDS_PC && rerun; outterrounds++ )
1959  {
1960  int nouterelims = 0;
1961  const SCIP_Bool beFast = (outterrounds == 0 && advancedrun);
1962  rerun = FALSE;
1963 
1964  SCIP_CALL( redLoopInnerPc(scip, g, redsollocal, dheap, vnoi, path, nodearrreal, heap, state,
1965  vbase, nodearrint, edgearrint, nodearrint2, nodearrchar, fixed, randnumgen, prizesum,
1966  dualascent, bred, reductbound, userec, nodereplacing, usestrongreds, beFast, &nouterelims) );
1967 
1968  if( advancedrun && g->terms > 2 )
1969  {
1970  const RPDA paramsda = {
1972  .useRec = userec,
1973  .extredMode = (nadvruns == 0 ? extred_fast : extred_full),
1974  .nodereplacing = nodereplacing,
1975  .useSlackPrune = FALSE,
1976  .pcmw_solbasedda = TRUE, .pcmw_useMultRoots = TRUE, .pcmw_markroots = TRUE, .pcmw_fastDa = FALSE };
1977 
1978  int danelims = 0;
1979  int implnelims = 0;
1980 
1981  degnelims = 0;
1982  advancedrun = FALSE;
1983  nadvruns++;
1984 
1985  if( rpc )
1986  {
1987  if( trySpgTrans && graph_pc_nNonLeafTerms(g) == 0 && graph_pc_nProperPotentialTerms(g) == 0 )
1988  {
1989  isSpg = TRUE;
1990  break;
1991  }
1992 
1993  SCIP_CALL( reduce_impliedProfitBasedRpc(scip, g, redsollocal, fixed, &implnelims) );
1994  if( verbose ) printf("implnelims: %d \n", implnelims);
1995 
1996  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, fixed, &danelims, randnumgen) );
1997  }
1998  else
1999  {
2000  reduce_sollocalSetOffset(*fixed, redsollocal);
2001  SCIP_CALL( reduce_sollocalRebuildTry(scip, g, redsollocal) );
2002  SCIP_CALL( reduce_daPcMw(scip, g, &paramsda, redsollocal, vnoi, nodearrreal, vbase, heap,
2003  state, nodearrchar, &danelims, randnumgen, prizesum) );
2004  }
2005 
2006  SCIP_CALL( reduce_simple_pc(scip, NULL, g, fixed, &degnelims, NULL, solnode) );
2007  nouterelims += danelims + degnelims + implnelims;
2008 
2009  if( nouterelims + implnelims > reductbound_global )
2010  {
2011  if( danelims + implnelims > reductbound_global )
2012  advancedrun = TRUE;
2013 
2014  rerun = TRUE;
2015  }
2016  }
2017 
2018  if( !rpc && tryrpc && g->terms > 2 )
2019  {
2020  assert(graph_pc_term2edgeIsConsistent(scip, g));
2021 
2022 #ifdef WITH_UG
2023  SCIP_CALL(graph_transPcmw2rooted(scip, g, prizesum, FALSE));
2024 #else
2025  SCIP_CALL(graph_transPcmw2rooted(scip, g, prizesum, TRUE));
2026 #endif
2027  rpc = (g->stp_type == STP_RPCSPG);
2028 
2029  if( rpc )
2030  {
2031  SCIP_CALL(reduce_unconnectedRpcRmw(scip, g, fixed));
2032  rerun = TRUE;
2033  advancedrun = dualascent;
2034  outterrounds = 1;
2035 
2036  if( trySpgTrans && graph_pc_nNonLeafTerms(g) == 0 && graph_pc_nProperPotentialTerms(g) == 0 )
2037  {
2038  isSpg = TRUE;
2039  break;
2040  }
2041  }
2042  }
2043  }
2044 
2045  if( dualascent && tryrpc)
2046  {
2047  //reduce_simple_aritculations(scip, g, NULL, NULL);
2049  }
2050 
2051  assert(graph_pc_term2edgeIsConsistent(scip, g));
2052  graph_pc_2trans(scip, g);
2053  graph_pc_presolExit(scip, g);
2054 
2055  graph_heap_free(scip, TRUE, TRUE, &dheap);
2056  SCIPfreeRandom(scip, &randnumgen);
2057 
2058  reduce_solFinalizeLocal(scip, g, redsol);
2059 
2060  if( isSpg )
2061  {
2062 #ifndef WITH_UG
2063  printf("changing problem to SPG \n");
2064  graph_printInfo(g);
2065 #endif
2067  }
2068 
2069  SCIPdebugMessage("Reduction Level PC 1: Fixed Cost = %.12e\n", *fixed);
2070  return SCIP_OKAY;
2071 }
2072 
2073 /** STP loop */
2075  SCIP* scip, /**< SCIP data structure */
2076  GRAPH* g, /**< graph data structure */
2077  REDBASE* redbase /**< parameters */
2078 )
2079 {
2080  REDSOLLOCAL* redsollocal;
2081  SCIP_RANDNUMGEN* randnumgen;
2082  const RPARAMS* redparameters = redbase->redparameters;
2083  SCIP_Real timelimit;
2084  const int reductbound = redparameters->reductbound;
2085  const SCIP_Bool fullreduce = redparameters->fullreduce;
2086  const SCIP_Bool nodereplacing = redparameters->nodereplacing;
2087  int dummy = 0;
2088  int nruns = 0;
2089  SCIP_Bool rerun = TRUE;
2090 #ifdef USE_FULLSEPA
2091  SCIP_Bool runFullSepa = TRUE;
2092 #endif
2093 
2094  assert(scip && g && redbase);
2095  assert(reductbound > 0);
2096  assert(graph_valid(scip, g));
2097  assert(graph_typeIsSpgLike(g));
2098 
2099  SCIP_CALL( reduce_solInitLocal(scip, g, redbase->redsol, &redsollocal));
2100  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
2101  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
2102 
2103  SCIP_CALL( reduce_contract0Edges(scip, g, redbaseGetSolnode(redsollocal, redbase), TRUE) );
2104  SCIP_CALL( reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &dummy, NULL) );
2105 
2106  /* reduction loop */
2107  do
2108  {
2109  SCIP_Bool wasDecomposed;
2110  rerun = FALSE;
2111 
2112  SCIP_CALL( redLoopInnerStp(scip, randnumgen, g, redsollocal, redbase, &wasDecomposed) );
2113 
2114  if( wasDecomposed )
2115  {
2116  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &dummy, NULL));
2117  break;
2118  }
2119 
2120  if( fullreduce && !SCIPisStopped(scip) )
2121  {
2122  int extendedelims = 0;
2123  const RPDA paramsda = {
2124  .damode = (nruns == 0 ? STP_DAMODE_FAST : STP_DAMODE_MEDIUM),
2125  .useRec = redparameters->userec,
2126  .extredMode = (nruns == 0 ? extred_fast : extred_full),
2127  .useSlackPrune = FALSE,
2128  .nodereplacing = nodereplacing};
2129 
2130  if( SCIPgetTotalTime(scip) > timelimit )
2131  break;
2132 
2133  if( redparameters->userec && (SCIPStpEnumRelaxIsPromising(g) || dpterms_isPromisingFully(g)) )
2134  break;
2135 
2136  assert(!rerun);
2137 
2138  SCIP_CALL( reduce_da(scip, g, &paramsda, redsollocal, redbaseGetOffsetPointer(redbase), &extendedelims, randnumgen) );
2139  reduceStatsPrint(fullreduce, "ext", extendedelims);
2140 
2141  SCIP_CALL(reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &extendedelims, NULL));
2142 
2143  if( nodereplacing )
2144  {
2145  int conflictnelims = 0;
2146 
2147  SCIP_CALL( reduce_fixedConflicts(scip, NULL, g, &conflictnelims) );
2148  reduceStatsPrint(fullreduce, "fixedconflict", conflictnelims);
2149  extendedelims += conflictnelims;
2150  }
2151 
2152  if( extendedelims > STP_RED_EXPENSIVEFACTOR * reductbound )
2153  {
2154  rerun = TRUE;
2155  }
2156 
2157 #ifdef USE_FULLSEPA
2158  // todo tune!
2159  if( !rerun && runFullSepa )
2160  {
2161  int sepanelims = 0;
2162 
2163  if( dpterms_isPromisingFully(g) )
2164  break;
2165 
2166  SCIP_CALL( reduce_termsepaFull(scip, g, redbaseGetSolnode(redsollocal, redbase), redbase, &sepanelims) );
2167 
2168  if( sepanelims > STP_RED_EXPENSIVEFACTOR * reductbound )
2169  SCIP_CALL( redLoopInnerStp(scip, randnumgen, g, redsollocal, redbase, &wasDecomposed) );
2170  else if( sepanelims > reductbound )
2171  SCIP_CALL( reduce_simple(scip, g, redbaseGetOffsetPointer(redbase), redbaseGetSolnode(redsollocal, redbase), &dummy, NULL));
2172 
2173  runFullSepa = FALSE;
2174  }
2175 #endif
2176  }
2177 
2178  nruns++;
2179  }
2180  while( rerun && !SCIPisStopped(scip) ); /* reduction loop */
2181 
2182  assert(graph_valid_ancestors(scip, g));
2183  SCIPfreeRandom(scip, &randnumgen);
2184 
2185  reduce_solFinalizeLocal(scip, g, redbase->redsol);
2186 
2187  return SCIP_OKAY;
2188 }
2189 
2190 
2191 /** reduces the graph */
2193  SCIP* scip, /**< SCIP data structure */
2194  GRAPH* graph, /**< graph structure */
2195  REDSOL* redsol, /**< primal solution container */
2196  int reductionlevel, /**< reduction level 0: none, 1: basic, 2: advanced */
2197  int minelims, /**< minimal amount of reductions to reiterate reduction methods */
2198  SCIP_Bool userec /**< use recombination heuristic? */
2199  )
2200 {
2201  int stp_type;
2202  SCIP_Real* offset = reduce_solGetOffsetPointer(redsol);
2203  SCIP_Bool graphWasInitialized = FALSE;
2204 
2205  assert(scip && offset);
2206  assert(graph_getFixedges(graph) == NULL);
2207  assert(reductionlevel == STP_REDUCTION_NONE || reductionlevel == STP_REDUCTION_BASIC || reductionlevel == STP_REDUCTION_ADVANCED );
2208  assert(minelims >= 0);
2209  assert(EQ(*offset, 0.0));
2210 
2211  stp_type = graph->stp_type;
2212 
2213  if( !graph_isSetUp(graph) )
2214  {
2215  SCIP_CALL( graph_initHistory(scip, graph) );
2216  }
2217 
2218  if( !graph_path_exists(graph) )
2219  {
2220  SCIP_CALL( graph_path_init(scip, graph) );
2221  graphWasInitialized = TRUE;
2222  }
2223 
2224  /* start with trivial preprocessing step */
2225  SCIP_CALL( reduce_unconnected(scip, graph) );
2226 
2227  /* if no reduction methods available for given problem, return */
2228  if( graph->stp_type == STP_BRMWCSP )
2229  {
2230  graph_path_exit(scip, graph);
2231  return SCIP_OKAY;
2232  }
2233 
2234  if( reductionlevel == STP_REDUCTION_BASIC )
2235  {
2236  if( stp_type == STP_PCSPG || stp_type == STP_RPCSPG )
2237  {
2238  SCIP_CALL( reduce_pc(scip, redsol, graph, minelims, FALSE, FALSE, TRUE, TRUE) );
2239  }
2240  else if( stp_type == STP_MWCSP || stp_type == STP_RMWCSP )
2241  {
2242  SCIP_CALL( reduce_mw(scip, redsol, graph, minelims, FALSE, FALSE, TRUE) );
2243  }
2244  else if( stp_type == STP_DHCSTP )
2245  {
2246  SCIP_CALL( reduce_hc(scip, graph, offset, minelims) );
2247  }
2248  else if( stp_type == STP_SAP || stp_type == STP_NWPTSPG )
2249  {
2250  SCIP_CALL( reduce_sap(scip, graph, FALSE, offset, minelims) );
2251  }
2252  else if( stp_type == STP_NWSPG )
2253  {
2254  SCIP_CALL( reduce_nw(scip, graph, offset, minelims) );
2255  }
2256  else if( stp_type == STP_DCSTP )
2257  {
2258  SCIP_CALL( reduce_dc(scip, graph, offset, minelims) );
2259  }
2260  else
2261  {
2262  assert(graph_typeIsSpgLike(graph));
2263 
2264  SCIP_CALL( reduce_stp(scip, graph, redsol, minelims, FALSE, TRUE, FALSE, TRUE) );
2265  }
2266  }
2267  else if( reductionlevel == STP_REDUCTION_ADVANCED )
2268  {
2269  if( stp_type == STP_PCSPG || stp_type == STP_RPCSPG )
2270  {
2271  SCIP_CALL( reduce_pc(scip, redsol, graph, minelims, TRUE, userec, TRUE, TRUE) );
2272 
2273  if( graph->stp_type == STP_SPG )
2274  {
2275  reduce_solReInitLocal(graph, redsol);
2276  SCIP_CALL( reduce_stp(scip, graph, redsol, minelims, TRUE, TRUE, userec, TRUE) );
2277  }
2278  }
2279  else if( stp_type == STP_MWCSP || stp_type == STP_RMWCSP )
2280  {
2281  SCIP_CALL( reduce_mw(scip, redsol, graph, minelims, TRUE, userec, TRUE) );
2282  }
2283  else if( stp_type == STP_DHCSTP )
2284  {
2285  SCIP_CALL( reduce_hc(scip, graph, offset, minelims) );
2286  }
2287  else if( stp_type == STP_SAP || stp_type == STP_NWPTSPG )
2288  {
2289  SCIP_CALL( reduce_sap(scip, graph, TRUE, offset, minelims) );
2290  }
2291  else if( stp_type == STP_NWSPG )
2292  {
2293  SCIP_CALL( reduce_nw(scip, graph, offset, minelims) );
2294  }
2295  else if( stp_type == STP_DCSTP )
2296  {
2297  SCIP_CALL( reduce_dc(scip, graph, offset, minelims) );
2298  }
2299  else
2300  {
2301  assert(graph_typeIsSpgLike(graph));
2302 
2303  SCIP_CALL( reduce_stp(scip, graph, redsol, minelims, TRUE, TRUE, userec, TRUE) );
2304  }
2305  }
2306 
2307  SCIP_CALL( reduce_unconnected(scip, graph) );
2308 
2309  /* NOTE: ugly, but necessary to allow for clean execution of decomposition routine during solving process */
2310  if( graph_typeIsSpgLike(graph) && userec )
2311  {
2312  int nartelims = 0;
2313  SCIP_CALL( reduce_articulations(scip, graph, offset, &nartelims) );
2314  }
2315 
2316  SCIPdebugMessage("offset : %f \n", *offset);
2317 
2318  assert(graph_valid(scip, graph));
2319 
2320  if( graphWasInitialized )
2321  graph_path_exit(scip, graph);
2322 
2323  return SCIP_OKAY;
2324 }
void SCIPfreeRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen)
SCIP_RETCODE graph_transPcmw2rooted(SCIP *, GRAPH *, SCIP_Real, SCIP_Bool)
Definition: graph_trans.c:809
int graph_get_nEdges(const GRAPH *)
Definition: graph_stats.c:548
SCIP_RETCODE reduce_ans(SCIP *, GRAPH *, int *)
Definition: reduce_alt.c:1470
int * edgearrint
Definition: reducedefs.h:114
SCIP_Real * nodearrreal
Definition: reducedefs.h:109
SCIP_RETCODE reduce_boundHopRc(SCIP *, GRAPH *, PATH *, SCIP_Real *, SCIP_Real *, SCIP_Real *, SCIP_Real, int *, int *, int *, int *, int *, SCIP_Bool)
Definition: reduce_bnd.c:1646
SCIP_Bool graph_valid_ancestors(SCIP *, const GRAPH *)
static volatile int nterms
Definition: interrupt.c:38
SCIP_RETCODE reduce_bdk(SCIP *, int, GRAPH *, int *)
SCIP_RETCODE reduce_da(SCIP *, GRAPH *, const RPDA *, REDSOLLOCAL *, SCIP_Real *, int *, SCIP_RANDNUMGEN *)
Definition: reduce_da.c:2471
SCIP_RETCODE reduce_simple_sap(SCIP *, GRAPH *, SCIP_Real *, int *)
SCIP_Bool graph_pc_isMw(const GRAPH *)
SCIP_RETCODE reduce_mw(SCIP *scip, REDSOL *redsol, GRAPH *g, int minelims, SCIP_Bool advanced, SCIP_Bool userec, SCIP_Bool usestrongreds)
Definition: reduce_base.c:1344
SCIP_RETCODE reduce_ansAdv2(SCIP *, GRAPH *, int *)
Definition: reduce_alt.c:1656
SCIP_RETCODE reduce_fixedConflicts(SCIP *, const int *, GRAPH *, int *)
SCIP_RETCODE reduce_nnp(SCIP *, GRAPH *, int *)
Definition: reduce_alt.c:2534
SCIP_Bool graph_pc_isRootedPcMw(const GRAPH *)
static SCIP_RETCODE execPc_NVSL(SCIP *scip, SCIP_Bool usestrongreds, GRAPH *g, PATH *vnoi, SCIP_Real *nodearrreal, SCIP_Real *fixed, int *edgearrint, int *vbase, int *neighb, int *distnode, int *solnode, STP_Bool *visited, int *nelims, int redbound, SCIP_Bool verbose, SCIP_Bool *rerun)
Definition: reduce_base.c:338
SCIP_RETCODE reduce_sdWalkTriangle(SCIP *, int, SCIP_Bool, GRAPH *, int *, SCIP_Real *, int *, STP_Bool *, DHEAP *, int *)
Definition: reduce_sd.c:3011
SCIP_Bool fullreduce
Definition: reducedefs.h:83
SCIP_RETCODE reduce_impliedProfitBased(SCIP *, int, GRAPH *, int *, SCIP_Real *, int *)
Definition: reduce_alt.c:2609
STP_REDTYPE
Definition: reduce_base.c:70
#define SCIPfreeMemoryArray(scip, ptr)
Definition: scip_mem.h:71
SCIP_Bool graph_valid(SCIP *, const GRAPH *)
Definition: graph_base.c:1480
SCIP_RETCODE SCIPgetRealParam(SCIP *scip, const char *name, SCIP_Real *value)
Definition: scip_param.c:298
SCIP_RETCODE reduce_redLoopStp(SCIP *scip, GRAPH *g, REDBASE *redbase)
Definition: reduce_base.c:2074
SCIP_RETCODE graph_pc_presolInit(SCIP *, GRAPH *)
Definition: graph_pcbase.c:815
void reduce_solFinalizeLocal(SCIP *, const GRAPH *, REDSOL *)
Definition: reduce_sol.c:735
int terms
Definition: graphdefs.h:192
static int * redbaseGetSolnode(REDSOLLOCAL *redsollocal, REDBASE *redbase)
Definition: reduce_base.c:417
void reduce_sollocalSetOffset(SCIP_Real, REDSOLLOCAL *)
Definition: reduce_sol.c:488
#define SCIPallocMemoryArray(scip, ptr, num)
Definition: scip_mem.h:55
int * nodearrint
Definition: reducedefs.h:113
int * nodearrint2
Definition: reducedefs.h:115
SCIP_Bool graph_pc_term2edgeIsConsistent(SCIP *, const GRAPH *)
Definition: graph_pcbase.c:874
SCIP_RETCODE reduce_nvAdv(SCIP *, const int *, GRAPH *, PATH *, SCIP_Real *, SCIP_Real *, int *, int *, int *, int *, int *)
Definition: reduce_alt.c:1121
void graph_pc_2org(SCIP *, GRAPH *)
SCIP_RETCODE graph_path_init(SCIP *, GRAPH *)
Definition: graph_path.c:480
#define FALSE
Definition: def.h:87
static SCIP_RETCODE execNvSl(SCIP *scip, const int *edgestate, GRAPH *g, PATH *vnoi, SCIP_Real *nodearrreal, SCIP_Real *fixed, int *edgearrint, int *vbase, int *neighb, int *distnode, int *solnode, STP_Bool *visited, int *nelims, int minelims)
Definition: reduce_base.c:196
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
#define STP_SAP
Definition: graphdefs.h:43
#define STP_DCSTP
Definition: graphdefs.h:47
Dynamic programming solver for Steiner tree (sub-) problems with small number of terminals.
static SCIP_RETCODE execPc_SD(SCIP *scip, GRAPH *g, PATH *vnoi, int *heap, int *state, int *vbase, int *nodesid, int *nodesorg, int *nelims, int redbound, SCIP_Bool verbose, SCIP_Bool *rerun)
Definition: reduce_base.c:274
SCIP_Bool reduce_solUsesNodesol(const REDSOL *)
Definition: reduce_sol.c:1236
SCIP_RETCODE reduce_sollocalInit(SCIP *, const GRAPH *, REDSOLLOCAL **)
Definition: reduce_sol.c:447
#define STP_SPG
Definition: graphdefs.h:42
#define SCIPdebugMessage
Definition: pub_message.h:87
#define FARAWAY
Definition: graphdefs.h:89
static SCIP_RETCODE redLoopInnerPc(SCIP *scip, GRAPH *g, REDSOLLOCAL *redsollocal, DHEAP *dheap, PATH *vnoi, PATH *path, SCIP_Real *nodearrreal, int *heap, int *state, int *vbase, int *nodearrint, int *edgearrint, int *nodearrint2, STP_Bool *nodearrchar, SCIP_Real *fixed, SCIP_RANDNUMGEN *randnumgen, SCIP_Real prizesum, SCIP_Bool dualascent, SCIP_Bool bred, int reductbound, SCIP_Bool userec, SCIP_Bool nodereplacing, SCIP_Bool usestrongreds, SCIP_Bool beFast, int *ninnerelims)
Definition: reduce_base.c:872
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:127
void graph_printInfo(const GRAPH *)
Definition: graph_stats.c:299
SCIP_RETCODE reduce_sollocalRebuildTry(SCIP *, GRAPH *, REDSOLLOCAL *)
Definition: reduce_sol.c:507
SCIP_RETCODE reduce_npv(SCIP *, GRAPH *, PATH *, int *, int *, int *, int *, int)
Definition: reduce_alt.c:2067
static SCIP_RETCODE redLoopInnerMw(SCIP *scip, GRAPH *g, REDSOLLOCAL *redsollocal, PATH *vnoi, SCIP_Real *nodearrreal, int *state, int *vbase, int *nodearrint, STP_Bool *nodearrchar, SCIP_Real *fixed, STP_Bool dualascent, STP_Bool bred, int redbound, SCIP_Bool userec, SCIP_Real prizesum)
Definition: reduce_base.c:691
SCIP_RETCODE reduce_solInitLocal(SCIP *, const GRAPH *, REDSOL *, REDSOLLOCAL **)
Definition: reduce_sol.c:717
SCIP_RETCODE reduce_daPcMw(SCIP *, GRAPH *, const RPDA *, REDSOLLOCAL *, PATH *, SCIP_Real *, int *, int *, int *, STP_Bool *, int *, SCIP_RANDNUMGEN *, SCIP_Real)
Definition: reduce_da.c:3174
#define STP_DAMODE_MEDIUM
Definition: reducedefs.h:45
int reduce_getMinNreductions(const GRAPH *g, int lowerbound)
Definition: reduce_base.c:1087
SCIP_RETCODE reduce_termsepaFull(SCIP *, GRAPH *, int *, REDBASE *, int *)
void reduce_baseFree(SCIP *scip, REDBASE **redbase)
Definition: reduce_base.c:1155
#define STP_REDUCTION_BASIC
Definition: reducedefs.h:40
SCIP_RETCODE reduce_sd(SCIP *, GRAPH *, REDBASE *, int *)
Definition: reduce_sd.c:1517
SCIP_RETCODE reduce_simple_hc(SCIP *, GRAPH *, SCIP_Real *, int *)
#define STP_REDUCTION_ADVANCED
Definition: reducedefs.h:41
#define STP_REDBOUND_PCBDK2
Definition: reduce_base.c:33
miscellaneous methods used for solving Steiner problems
#define STP_REDBOUND_PCBDK
Definition: reduce_base.c:32
void graph_pc_2trans(SCIP *, GRAPH *)
int graph_get_nTerms(const GRAPH *)
Definition: graph_stats.c:560
SCIP_RETCODE reduce_chain2(SCIP *, GRAPH *, PATH *, int *, int *, int *, int *, int)
Definition: reduce_alt.c:2445
SCIP_RETCODE graph_transRpc2SpgTrivial(SCIP *, GRAPH *)
Definition: graph_trans.c:505
SCIP_RETCODE reduce_redLoopPc(SCIP *scip, REDSOL *redsol, GRAPH *g, PATH *vnoi, PATH *path, SCIP_Real *nodearrreal, int *heap, int *state, int *vbase, int *nodearrint, int *edgearrint, int *nodearrint2, STP_Bool *nodearrchar, SCIP_Bool dualascent, SCIP_Bool bred, SCIP_Bool tryrpc, int reductbound, SCIP_Bool userec, SCIP_Bool nodereplacing, SCIP_Bool usestrongreds)
Definition: reduce_base.c:1896
int stp_type
Definition: graphdefs.h:264
SCIP_RETCODE reduce_simple_pc(SCIP *, const int *, GRAPH *, SCIP_Real *, int *, int *, int *)
void graph_get_nVET(const GRAPH *, int *, int *, int *)
Definition: graph_stats.c:571
int * reduce_sollocalGetSolnode(REDSOLLOCAL *)
Definition: reduce_sol.c:651
SCIP_Bool usestrongreds
Definition: reducedefs.h:84
static SCIP_RETCODE execPc_BDk(SCIP *scip, GRAPH *g, PATH *pathtail, PATH *pathhead, int *heap, int *statetail, int *statehead, int *memlbltail, int *memlblhead, int *nelims, int limit, SCIP_Real *offset, int redbound, SCIP_Bool verbose, SCIP_Bool *rerun)
Definition: reduce_base.c:307
SCIP_Bool boundreduce
Definition: reducedefs.h:78
#define SCIPfreeBufferArrayNull(scip, ptr)
Definition: scip_mem.h:128
SCIP_RETCODE reduce_nw(SCIP *scip, GRAPH *g, SCIP_Real *fixed, int minelims)
Definition: reduce_base.c:1669
SCIP_RETCODE graph_initHistory(SCIP *, GRAPH *)
Definition: graph_base.c:695
SCIP_RETCODE reduce_sl(SCIP *, const int *, GRAPH *, PATH *, SCIP_Real *, int *, int *, STP_Bool *, int *, int *)
Definition: reduce_alt.c:665
void print(const Container &container, const std::string &prefix="", const std::string &suffix="", std::ostream &os=std::cout, bool negate=false, int prec=6)
SCIP_RETCODE reduce_unconnected(SCIP *, GRAPH *)
SCIP_Bool dualascent
Definition: reducedefs.h:77
#define NULL
Definition: lpi_spx1.cpp:155
PC_REDTYPE
Definition: reduce_base.c:69
SCIP_RETCODE reduce_sdspSap(SCIP *, GRAPH *, PATH *, PATH *, int *, int *, int *, int *, int *, int *, int)
Definition: reduce_sd.c:2661
Steiner tree enumeration relaxator.
#define STP_DHCSTP
Definition: graphdefs.h:52
SCIP_Real reduce_solGetOffset(const REDSOL *)
Definition: reduce_sol.c:1137
#define STP_RMWCSP
Definition: graphdefs.h:54
#define EQ(a, b)
Definition: portab.h:79
SCIP_RETCODE reduce_bound(SCIP *, GRAPH *, PATH *, SCIP_Real *, SCIP_Real *, SCIP_Real *, int *, int *, int *, int *)
Definition: reduce_bnd.c:655
int knots
Definition: graphdefs.h:190
SCIP_RETCODE reduce_baseInit(SCIP *scip, const GRAPH *g, REDBASE **redbase)
Definition: reduce_base.c:1120
#define SCIP_CALL(x)
Definition: def.h:384
SCIP_RETCODE reduce_sap(SCIP *scip, GRAPH *g, SCIP_Bool dualascent, SCIP_Real *fixed, int minelims)
Definition: reduce_base.c:1557
#define STP_PCSPG
Definition: graphdefs.h:44
#define STP_REDBOUND_SDSP2
Definition: reduce_base.c:31
#define STP_RED_MAXNINNROUNDS
Definition: reduce_base.c:44
SCIP_RETCODE reduce_stp(SCIP *scip, GRAPH *g, REDSOL *redsol, int minelims, SCIP_Bool dualascent, SCIP_Bool nodereplacing, SCIP_Bool userec, SCIP_Bool usestrongreds)
Definition: reduce_base.c:1181
unsigned char STP_Bool
Definition: portab.h:34
#define STP_DAMODE_FAST
Definition: reducedefs.h:44
static void reduceStatsPrint(SCIP_Bool print, const char *method, int nelims)
Definition: reduce_base.c:178
SCIP_RETCODE SCIPcreateRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen, unsigned int initialseed, SCIP_Bool useglobalseed)
SCIP_RETCODE reduce_hc(SCIP *scip, GRAPH *g, SCIP_Real *fixed, int minelims)
Definition: reduce_base.c:1403
SCIP_RETCODE reduce_boundHop(SCIP *, GRAPH *, PATH *, SCIP_Real *, SCIP_Real *, SCIP_Real *, int *, int *, int *, int *)
Definition: reduce_bnd.c:1323
SCIP_RETCODE reduce_pc(SCIP *scip, REDSOL *redsol, GRAPH *g, int minelims, SCIP_Bool advanced, SCIP_Bool userec, SCIP_Bool nodereplacing, SCIP_Bool usestrongreds)
Definition: reduce_base.c:1260
SCIP_RETCODE reduce_simple(SCIP *, GRAPH *, SCIP_Real *, int *, int *, int *)
propagator for Steiner tree problems, using the LP reduced costs
SCIP_RETCODE reduce_bd34(SCIP *, GRAPH *, PATH *, PATH *, int *, int *, int *, int *, int *, int *, int, SCIP_Real *)
Definition: reduce_sd.c:4516
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:115
void graph_heap_free(SCIP *, SCIP_Bool, SCIP_Bool, DHEAP **)
Definition: graph_util.c:1034
SCIP_RETCODE reduce_sdPc(SCIP *, GRAPH *, PATH *, int *, int *, int *, int *, int *, int *)
Definition: reduce_sd.c:2023
#define STP_RED_GLBFACTOR
Definition: reduce_base.c:41
#define SCIP_Bool
Definition: def.h:84
SCIP_RETCODE reduce_deleteConflictEdges(SCIP *, GRAPH *)
Definition: reduce_ext.c:791
#define STP_NWPTSPG
Definition: graphdefs.h:48
static int getWorkLimitsPc(const GRAPH *g, int roundnumber, SCIP_Bool beFast, enum PC_REDTYPE redtype)
Definition: reduce_base.c:118
#define STP_REDBOUND_BDK
Definition: reduce_base.c:35
static SCIP_Real * redbaseGetOffsetPointer(REDBASE *redbase)
Definition: reduce_base.c:403
static int getWorkLimitsStp(const GRAPH *g, int roundnumber, SCIP_Bool fullreduce, enum STP_REDTYPE redtype)
Definition: reduce_base.c:75
SCIP_RETCODE reduce_rpt(SCIP *, GRAPH *, SCIP_Real *, int *)
SCIP_RETCODE reduce_pathreplace(SCIP *, GRAPH *, int *)
Definition: reduce_path.c:1119
STP_Bool * nodearrchar
Definition: reducedefs.h:116
REDSOL * redsol
Definition: reducedefs.h:105
#define MAX(x, y)
Definition: tclique_def.h:83
SCIP_RETCODE reduce_ansAdv(SCIP *, GRAPH *, int *, SCIP_Bool)
Definition: reduce_alt.c:1544
SCIP_RETCODE reduce_dapaths(SCIP *, GRAPH *, SCIP_Real *, int *)
Definition: reduce_da.c:2397
void reduce_solReInitLocal(const GRAPH *, REDSOL *)
Definition: reduce_sol.c:793
int graph_pc_nProperPotentialTerms(const GRAPH *)
#define STP_BND_THRESHOLD
Definition: reduce_base.c:48
IDX * graph_getFixedges(const GRAPH *)
SCIP_RETCODE reduce_sdsp(SCIP *, GRAPH *, PATH *, int *, int *, int *, int *, int *, int *, int, SCIP_Bool)
Definition: reduce_sd.c:4020
SCIP_RETCODE reduce_unconnectedForDirected(SCIP *, GRAPH *)
SCIP_RETCODE reduce_boundHopR(SCIP *, GRAPH *, PATH *, SCIP_Real *, SCIP_Real *, SCIP_Real *, int *, int *, int *, int *, int *)
Definition: reduce_bnd.c:1518
void graph_path_exit(SCIP *, GRAPH *)
Definition: graph_path.c:515
SCIP_Bool graph_pc_isPc(const GRAPH *)
SCIP_Bool dpterms_isPromisingFully(const GRAPH *)
Definition: dpterms_base.c:563
Portable definitions.
int graph_pc_nNonLeafTerms(const GRAPH *)
SCIP_Bool graph_typeIsUndirected(const GRAPH *)
Definition: graph_stats.c:69
#define SCIPfreeMemory(scip, ptr)
Definition: scip_mem.h:69
#define STP_NWSPG
Definition: graphdefs.h:46
static SCIP_RETCODE execPc_BND(SCIP *scip, GRAPH *graph, PATH *vnoi, SCIP_Real *radius, SCIP_Real *offset, int *heap, int *state, int *vbase, int *nelims, int redbound, SCIP_Bool verbose, SCIP_Bool *rerun)
Definition: reduce_base.c:371
void reduce_sollocalFree(SCIP *, REDSOLLOCAL **)
Definition: reduce_sol.c:472
#define STP_RED_EDGELIMIT_HUGE
Definition: reduce_base.c:43
#define STP_REDBOUND_BDK2
Definition: reduce_base.c:36
SCIP_Bool SCIPisGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_RETCODE reduce_unconnectedRpcRmw(SCIP *, GRAPH *, SCIP_Real *)
SCIP_Real * cost
Definition: graphdefs.h:221
#define STP_REDUCTION_NONE
Definition: reducedefs.h:39
#define STP_REDBOUND_SDSTAR2
Definition: reduce_base.c:38
SCIP_RETCODE graph_heap_create(SCIP *, int, int *, DENTRY *, DHEAP **)
Definition: graph_util.c:992
SCIP_RETCODE reduce_simple_mw(SCIP *, GRAPH *, int *, SCIP_Real *, int *)
SCIP_Real * reduce_solGetOffsetPointer(REDSOL *)
Definition: reduce_sol.c:1285
#define SCIP_Real
Definition: def.h:177
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:694
SCIP_RETCODE reduce_sdStarBiased(SCIP *, int, SCIP_Bool, GRAPH *, int *)
Definition: reduce_sd.c:3332
SCIP_RETCODE reduce_sdStarPc2(SCIP *, int, SCIP_Bool, GRAPH *, SCIP_Real *, int *, int *, STP_Bool *, DHEAP *, int *)
Definition: reduce_sd.c:3394
#define STP_BRMWCSP
Definition: graphdefs.h:55
shortest paths based primal heuristics for Steiner problems
SCIP_RETCODE reduce_sdWalkExt(SCIP *, int, SCIP_Bool, GRAPH *, SCIP_Real *, int *, int *, int *, STP_Bool *, int *)
Definition: reduce_sd.c:3827
int edges
Definition: graphdefs.h:219
SCIP_Bool SCIPStpEnumRelaxIsPromising(const GRAPH *graph)
SCIP_Real SCIPgetTotalTime(SCIP *scip)
Definition: scip_timing.c:342
int graph_get_nNodes(const GRAPH *)
Definition: graph_stats.c:536
SCIP_RETCODE reduce_boundHopDa(SCIP *, GRAPH *, int *, SCIP_RANDNUMGEN *)
Definition: reduce_bnd.c:1286
BIDECPARAMS * bidecompparams
Definition: reducedefs.h:103
#define STP_RED_EXPENSIVEFACTOR
Definition: reduce_base.c:40
#define SCIPallocMemory(scip, ptr)
Definition: scip_mem.h:51
SCIP_Bool graph_isSetUp(const GRAPH *)
Definition: graph_base.c:1240
static SCIP_RETCODE redLoopInnerStp(SCIP *scip, SCIP_RANDNUMGEN *randnumgen, GRAPH *g, REDSOLLOCAL *redsollocal, REDBASE *redbase, SCIP_Bool *wasDecomposed)
Definition: reduce_base.c:441
RPARAMS * redparameters
Definition: reducedefs.h:102
#define STP_RED_MWTERMBOUND
Definition: reduce_base.c:39
SCIP_Bool SCIPisLE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_RETCODE reduce_dc(SCIP *scip, GRAPH *g, SCIP_Real *fixed, int minelims)
Definition: reduce_base.c:1729
#define nnodes
Definition: gastrans.c:65
SCIP_RETCODE reduce_bidecomposition(SCIP *, GRAPH *, REDBASE *, int *, SCIP_Bool *)
Definition: reduce_sepa.c:1039
SCIP_RETCODE reduce_impliedProfitBasedRpc(SCIP *, GRAPH *, REDSOLLOCAL *, SCIP_Real *, int *)
Definition: reduce_alt.c:2664
#define STP_REDBOUND_SDSTAR
Definition: reduce_base.c:37
SCIP_RETCODE reduce_boundMw(SCIP *, GRAPH *, PATH *, SCIP_Real *, int *, int *, int *, int *, int *)
Definition: reduce_bnd.c:1058
#define STP_REDBOUND_SDSP
Definition: reduce_base.c:30
#define STP_RED_MAXNOUTROUNDS_PC
Definition: reduce_base.c:45
SCIP_RETCODE reduce_contract0Edges(SCIP *, GRAPH *, int *, SCIP_Bool)
SCIP_RETCODE reduce_articulations(SCIP *, GRAPH *, SCIP_Real *, int *)
Definition: reduce_sepa.c:1130
Minimum cut routines for Steiner problems.
SCIP_Real graph_pc_getPosPrizeSum(SCIP *, const GRAPH *)
SCIP_Bool graph_typeIsSpgLike(const GRAPH *)
Definition: graph_stats.c:57
SCIP_Bool nodereplacing
Definition: reducedefs.h:79
void reduce_simple_dc(SCIP *, GRAPH *)
SCIP_RETCODE reduce_exec(SCIP *scip, GRAPH *graph, REDSOL *redsol, int reductionlevel, int minelims, SCIP_Bool userec)
Definition: reduce_base.c:2192
includes various reduction methods for Steiner tree problems
SCIP_RETCODE reduce_redLoopMw(SCIP *scip, REDSOL *redsol, GRAPH *g, PATH *vnoi, SCIP_Real *nodearrreal, int *state, int *vbase, int *nodearrint, STP_Bool *nodearrchar, STP_Bool advanced, STP_Bool bred, STP_Bool tryrmw, int redbound, SCIP_Bool userec, SCIP_Bool usestrongreds)
Definition: reduce_base.c:1771
SCIP_Real reduce_sollocalGetUpperBound(const REDSOLLOCAL *)
Definition: reduce_sol.c:633
#define STP_RPCSPG
Definition: graphdefs.h:45
#define STP_MWCSP
Definition: graphdefs.h:51
#define STP_RED_EXTENSIVE
Definition: reduce_base.c:29
void graph_pc_presolExit(SCIP *, GRAPH *)
Definition: graph_pcbase.c:858
SCIP_Bool graph_path_exists(const GRAPH *)
Definition: graph_path.c:497
SCIP callable library.
SCIP_Bool reduce_sollocalUsesNodesol(const REDSOLLOCAL *)
Definition: reduce_sol.c:554