dec_stairheur.cpp
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program */
4 /* GCG --- Generic Column Generation */
5 /* a Dantzig-Wolfe decomposition based extension */
6 /* of the branch-cut-and-price framework */
7 /* SCIP --- Solving Constraint Integer Programs */
8 /* */
9 /* Copyright (C) 2010-2018 Operations Research, RWTH Aachen University */
10 /* Zuse Institute Berlin (ZIB) */
11 /* */
12 /* This program is free software; you can redistribute it and/or */
13 /* modify it under the terms of the GNU Lesser General Public License */
14 /* as published by the Free Software Foundation; either version 3 */
15 /* of the License, or (at your option) any later version. */
16 /* */
17 /* This program is distributed in the hope that it will be useful, */
18 /* but WITHOUT ANY WARRANTY; without even the implied warranty of */
19 /* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the */
20 /* GNU Lesser General Public License for more details. */
21 /* */
22 /* You should have received a copy of the GNU Lesser General Public License */
23 /* along with this program; if not, write to the Free Software */
24 /* Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.*/
25 /* */
26 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
27 
39 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
40 
41 #include <cassert>
42 #include <cstring>
43 #include <algorithm>
44 #include <vector>
45 #include <iostream>
46 
47 #include "dec_stairheur.h"
48 #include "cons_decomp.h"
49 #include "struct_decomp.h"
50 #include "pub_decomp.h"
51 #include "scip_misc.h"
52 #include "scip/pub_misc.h"
53 #include "scip/struct_var.h"
54 #include "gcg.h"
55 #include "class_seeed.h"
56 #include "class_seeedpool.h"
57 #include "scip/clock.h"
58 
59 #define DEC_DETECTORNAME "stairheur"
60 #define DEC_DESC "detects staircase matrices via matrix reordering"
61 #define DEC_PRIORITY 1200
62 #define DEC_FREQCALLROUND 1
63 #define DEC_MAXCALLROUND INT_MAX
64 #define DEC_MINCALLROUND 0
65 #define DEC_FREQCALLROUNDORIGINAL 1
66 #define DEC_MAXCALLROUNDORIGINAL INT_MAX
67 #define DEC_MINCALLROUNDORIGINAL 0
68 #define DEC_DECCHAR 's'
69 #define DEC_ENABLED FALSE
70 #define DEC_ENABLEDORIGINAL FALSE
71 #define DEC_ENABLEDFINISHING FALSE
72 #define DEC_ENABLEDPOSTPROCESSING FALSE
73 #define DEC_SKIP FALSE
74 #define DEC_USEFULRECALL FALSE
75 #define DEC_LEGACYMODE FALSE
77 /* Default parameter settings*/
78 #define DEFAULT_NCONSSPERBLOCK 32
79 #define DEFAULT_MAXBLOCKS 20
80 #define DEFAULT_MINBLOCKS 2
81 #define DEFAULT_DESIREDBLOCKS 0
82 #define DEFAULT_DYNAMICBLOCKING FALSE
83 #define DEFAULT_STATICBLOCKING TRUE
84 #define DEFAULT_BLOCKINGASSOONASPOSSIBLE FALSE
85 #define DEFAULT_MULTIPLEDECOMPS TRUE
86 #define DEFAULT_MAXITERATIONSROC 1000000
88 using std::find;
89 using std::vector;
90 using std::swap;
91 
92 //#define WRITEALLOUTPUT
93 
94 
95 #ifdef WRITEALLOUTPUT
96 #define DWSOLVER_REFNAME(name, blocks, cons, dummy) "%s_%d_%d_%.1f_ref.txt", (name), (blocks), (cons), (dummy)
97 #define GP_NAME(name, blocks, cons, dummy) "%s_%d_%d_%.1f_%d.gp", (name), (blocks), (cons), (dummy)
98 #endif
99 /*
100  * Data structures
101  */
102 
104 struct IndexMap
105 {
106  SCIP_HASHMAP* indexcons;
107  SCIP_HASHMAP* consindex;
108  SCIP_HASHMAP* indexvar;
109  SCIP_HASHMAP* varindex;
110 };
111 typedef struct IndexMap INDEXMAP;
112 
114 struct DEC_DetectorData
115 {
116  SCIP_HASHMAP* constoblock;
117  int blocks;
119  int maxblocks;
120  int minblocks;
122  int* ibegin;
123  int* iend;
124  int* jbegin;
125  int* jend;
126  int* jmin;
127  int* jmax;
128  int* minV;
129  int* width;
131  vector<int>* rowsWithConstrictions;
132  vector<int>* blockedAfterrow;
134  SCIP_Bool dynamicblocking;
135  SCIP_Bool staticblocking;
137  SCIP_Bool multipledecomps;
139 };
140 
141 /*
142  * Local methods
143  */
144 
145 
146 /* debugging methods */
148  vector<vector<int> > l
149  )
150 {
151  vector<int>::iterator inner;
152  vector<vector<int> >::iterator outer;
153  SCIPdebugPrintf("S:");
154  for( outer = l.begin(); outer != l.end(); ++outer)
155  {
156  SCIPdebugPrintf("\t");
157  for( inner = outer->begin(); inner != outer->end(); ++inner)
158  {
159  SCIPdebugPrintf(" %d", *inner);
160  }
161  SCIPdebugPrintf(".\n");
162  }
163  SCIPdebugPrintf("Done\n");
164 }
165 
167  vector<int > l
168  )
169 {
170  vector<int>::iterator inner;
171  for( inner = l.begin(); inner != l.end(); ++inner)
172  {
173  SCIPdebugPrintf(" %d", *inner);
174  }
175 }
176 
177 
188 static
190  int from,
191  int to
192  )
193 {
194  vector<int> v;
195 
196  for( int i = from; i <= to; ++i )
197  {
198  v.push_back(i);
199  }
200  return v;
201 }
202 
203 
210 static
212  vector<vector<int> > &l,
213  vector<int> order)
214 {
215 
216  vector<vector<int> > new_vector;
217  vector<int>::iterator it1;
218 
219  for( it1 = order.begin(); it1 != order.end(); ++it1 )
220  {
221  new_vector.push_back(l[(size_t)*it1-1]);
222 
223  }
224  l.swap(new_vector);
225 
226 }
227 
229 static
230 SCIP_RETCODE indexmapCreate(
231  SCIP* scip,
232  INDEXMAP** indexmap,
233  int nconss,
234  int nvars
235 )
236 {
237  INDEXMAP* imap = NULL;
238  assert(scip != NULL);
239  assert(nconss > 0);
240  assert(nvars > 0);
241  SCIP_CALL( SCIPallocMemory(scip, &imap) );
242  assert(imap != NULL);
243 
244  SCIP_CALL( SCIPhashmapCreate(&imap->indexvar, SCIPblkmem(scip), nvars) );
245  SCIP_CALL( SCIPhashmapCreate(&imap->varindex, SCIPblkmem(scip), nvars) );
246  SCIP_CALL( SCIPhashmapCreate(&imap->indexcons, SCIPblkmem(scip), nconss) );
247  SCIP_CALL( SCIPhashmapCreate(&imap->consindex, SCIPblkmem(scip), nconss) );
248 
249  *indexmap = imap;
250  return SCIP_OKAY;
251 }
252 
254 static
256  SCIP* scip,
257  INDEXMAP** indexmap
258  )
259 {
260  SCIPhashmapFree(&(*indexmap)->indexvar);
261  SCIPhashmapFree(&(*indexmap)->varindex);
262  SCIPhashmapFree(&(*indexmap)->indexcons);
263  SCIPhashmapFree(&(*indexmap)->consindex);
264  SCIPfreeMemory(scip, indexmap);
265 }
266 
268 static
269 SCIP_RETCODE indexmapInit(
270  INDEXMAP* indexmap,
271  SCIP_VAR** vars,
272  int nvars,
273  SCIP_CONS** conss,
274  int nconss,
275  int* hashmapindices
276  )
277 {
278  int i;
279  int* hashmapindex;
280 
281  for( i = 0; i < nvars; ++i )
282  {
283  SCIP_VAR* var = vars[i];
284  /* careful: hashmapindex+1, because '0' is treated as an empty hashmap entry, which causes an error */
285  hashmapindex = hashmapindices + i+1;
286  assert( ! SCIPhashmapExists(indexmap->indexvar, (void*) hashmapindex));
287  SCIP_CALL( SCIPhashmapInsert(indexmap->indexvar, (void*) hashmapindex, (void*) var) );
288  assert( ! SCIPhashmapExists(indexmap->varindex, (void*) var));
289  SCIP_CALL( SCIPhashmapInsert(indexmap->varindex, (void*) var, (void*) hashmapindex) );
290  }
291 
292  for( i = 0; i < nconss; ++i )
293  {
294  SCIP_CONS* cons = conss[i];
295  /* careful: hashmapindex+1, because '0' is treated as an empty hashmap entry, which causes an error */
296  hashmapindex = hashmapindices + i+1;
297  assert( ! SCIPhashmapExists(indexmap->indexcons, (void*) hashmapindex));
298  SCIP_CALL( SCIPhashmapInsert(indexmap->indexcons, (void*) hashmapindex, (void*) cons) );
299  assert( ! SCIPhashmapExists(indexmap->consindex, (void*) cons));
300  SCIP_CALL( SCIPhashmapInsert(indexmap->consindex, (void*) cons, (void*) hashmapindex) );
301  }
302 
303  return SCIP_OKAY;
304 }
305 
307 static
308 SCIP_RETCODE indexmapInit(
309  INDEXMAP* indexmap,
310  gcg::Seeed* seeed,
311  int* hashmapindices
312  )
313 {
314  int i;
315  int* hashmapindex;
316  int nvars = seeed->getNOpenvars();
317  int nconss = seeed->getNOpenconss();
318 
319  for( i = 0; i < nvars; ++i )
320  {
321  int var = seeed->getOpenvars()[i];
322  /* careful: hashmapindex+1, because '0' is treated as an empty hashmap entry, which causes an error */
323  hashmapindex = hashmapindices + i+1;
324  assert( ! SCIPhashmapExists(indexmap->indexvar, (void*) hashmapindex));
325  SCIP_CALL( SCIPhashmapInsert(indexmap->indexvar, (void*) hashmapindex, (void*)(size_t) var) );
326  assert( ! SCIPhashmapExists(indexmap->varindex, (void*)(size_t) var));
327  SCIP_CALL( SCIPhashmapInsert(indexmap->varindex, (void*)(size_t) var, (void*) hashmapindex) );
328  }
329 
330  for( i = 0; i < nconss; ++i )
331  {
332  int cons = seeed->getOpenconss()[i];
333  /* careful: hashmapindex+1, because '0' is treated as an empty hashmap entry, which causes an error */
334  hashmapindex = hashmapindices + i+1;
335  assert( ! SCIPhashmapExists(indexmap->indexcons, (void*) hashmapindex));
336  SCIP_CALL( SCIPhashmapInsert(indexmap->indexcons, (void*) hashmapindex, (void*)(size_t) cons) );
337  assert( ! SCIPhashmapExists(indexmap->consindex, (void*)(size_t) cons));
338  SCIP_CALL( SCIPhashmapInsert(indexmap->consindex, (void*)(size_t) cons, (void*) hashmapindex) );
339  }
340 
341  return SCIP_OKAY;
342 }
343 
344 /* debug ? */
345 #ifdef WRITEALLOUTPUT
346 
347 static const char* getProbNameWithoutPath(
348  SCIP* scip
349 )
350 {
351  const char* pname;
352  /* remove '/' from problem name */
353  pname = strrchr(SCIPgetProbName(scip), '/');
354  if( pname == NULL )
355  {
356  pname = SCIPgetProbName(scip);
357  }
358  else
359  {
360  pname = pname+1;
361  }
362  return pname;
363 }
364 
365 
366 static void checkConsistencyOfIndexarrays(DEC_DETECTORDATA* detectordata, int nvars, int nconss)
367 {
368  int i;
369  for( i = 0; i < nconss - 1; ++i )
370  {
371  assert(detectordata->ibegin[i] <= detectordata->ibegin[i+1]);
372  }
373  for( i = 0; i < nvars - 1; ++i )
374  {
375  assert(detectordata->jbegin[i] <= detectordata->jbegin[i+1]);
376  }
377 }
378 
379 
380 /* debug ? */
385 static
386 SCIP_RETCODE plotInitialProblem(
387  SCIP* scip,
388  DEC_DETECTORDATA* detectordata,
389  char* filename
390 )
391 {
392  FILE* output;
393  char datafile[256];
394  char gpfile[256];
395  char pdffile[256];
396 
397  int nconss;
398  nconss = SCIPgetNConss(scip);
399 
400  /* filenames */
401  sprintf(datafile, "%s.dat", filename);
402  sprintf(gpfile, "%s.gp", filename);
403  sprintf(pdffile, "%s.pdf", filename);
404  output = fopen(datafile, "w");
405  if (output == NULL)
406  {
407  SCIPinfoMessage(scip, NULL, "Can't open file for output in plotProblem!\n");
408  }
409  else
410  {
411  int i;
412 
413  for( i = 0; i < nconss; ++i )
414  {
415  int j;
416  SCIP_Bool success;
417  SCIP_VAR** curvars;
418  int ncurvars;
419  int consindex;
420  SCIP_CONS* cons = SCIPgetConss(scip)[i];
421  consindex = *((int*) SCIPhashmapGetImage(detectordata->indexmap->consindex, cons));
422  assert(SCIPhashmapExists(detectordata->indexmap->consindex, cons));
423 
424  /* Get array of variables from constraint */
425  SCIP_CALL( SCIPgetConsNVars(scip, cons, &ncurvars, &success) );
426  assert(success);
427  SCIP_CALL( SCIPallocMemoryArray(scip, &curvars, ncurvars) );
428  SCIP_CALL( SCIPgetConsVars(scip, cons, curvars, ncurvars, &success) );
429  assert(success);
430  for( j = 0; j < ncurvars; ++j )
431  {
432  SCIP_VAR* var = SCIPvarGetProbvar(curvars[j]);
433  assert(SCIPhashmapExists(detectordata->indexmap->varindex, var));
434  int varindex = *((int*) SCIPhashmapGetImage(detectordata->indexmap->varindex, var));
435  assert(varindex <= SCIPgetNVars(scip));
436  assert(varindex > 0);
437  assert(consindex <= SCIPgetNConss(scip));
438  assert(consindex > 0);
439  fprintf(output, "%d %d\n", varindex, consindex);
440  }
441  SCIPfreeMemoryArray(scip, &curvars);
442  }
443  }
444  fclose(output);
445 
446  /* write Gnuplot file */
447  output = fopen(gpfile, "w");
448  fprintf(output, "set terminal pdf\nset output \"%s\"\nunset xtics\nunset ytics\nunset border\nset pointsize 0.05\nset xrange [0:%i]\nset yrange[%i:0]\nplot '%s' lt 0 pt 5 notitle", pdffile, SCIPgetNVars(scip), nconss, datafile);
449  fclose(output);
450  return SCIP_OKAY;
451 }
452 
456 static
457 void plotMinV(
458  SCIP* scip,
459  DEC_DETECTORDATA* detectordata,
460  char* filename
461 )
462 {
463  FILE* output;
464  char datafile[256];
465  char blockingfile[256];
466  char gpfile[256];
467  char pdffile[256];
468  int nconss;
469  vector<int>::iterator it1;
470 
471  nconss = SCIPgetNConss(scip);
472 
473  /* filenames */
474  sprintf(datafile, "%s.dat", filename);
475  sprintf(blockingfile, "%s_blocked_at.dat", filename);
476  sprintf(gpfile, "%s.gp", filename);
477  sprintf(pdffile, "%s.pdf", filename);
478 
479  /* datafile */
480  output = fopen(datafile, "w");
481  if (output == NULL)
482  {
483  SCIPinfoMessage(scip, NULL, "Can't open file for output in plotMinV!\n");
484  }
485  else
486  {
487  int i;
488 
489  /* write data to datafile */
490  for( i = 0; i < nconss -1; ++i )
491  {
492  fprintf(output, "%d\n", detectordata->minV[i]);
493  }
494  }
495  fclose(output);
496 
497  /* blocking points */
498  output = fopen(blockingfile, "w");
499  if (output == NULL)
500  {
501  SCIPinfoMessage(scip, NULL, "Can't open file for blocking output in plotMinV!\n");
502  }
503  else
504  {
505  /* write data to blockingfile */
506  for( it1 = detectordata->blockedAfterrow->begin(); it1 != detectordata->blockedAfterrow->end(); ++it1 )
507  {
508  fprintf(output, "%d %d\n", *it1-1, detectordata->minV[*it1-1]);
509  }
510  }
511  fclose(output);
512  /* write Gnuplot file */
513  output = fopen(gpfile, "w");
514  fprintf(output, "set terminal pdf\nset output \"%s\"\nset style line 1 lt 1 lc rgb \"black\"\nplot '%s' title '# verb. Variablen' ls 1 with lines, \\\n '%s' lt 0 pt 4 with points title \"Blockgrenze\"", pdffile, datafile, blockingfile);
515  fclose(output);
516 }
517 
518 
519 
520 #endif
521 
522 
523 
525 static
526 SCIP_RETCODE initData(
527  SCIP* scip,
528  DEC_DETECTORDATA* detectordata
529  )
530 {
531  int i;
532  int nvars;
533  int nconss;
534 
535  assert(scip != NULL);
536  assert(detectordata != NULL);
537 
538  nvars = SCIPgetNVars(scip);
539  nconss = SCIPgetNConss(scip);
540  detectordata->maxblocks = MIN(nconss, detectordata->maxblocks);
541 
542  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->ibegin, nconss) );
543  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->iend, nconss) );
544  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jbegin, nvars) );
545  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jend, nvars) );
546  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jmin, nconss) );
547  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jmax, nconss) );
548  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->minV, nconss-1) );
549  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->width, nconss) );
550  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->hashmapindices, (size_t)MAX(nvars, nconss) + 1) );
551  for( i = 0; i < MAX(nvars, nconss)+1; ++i )
552  {
553  detectordata->hashmapindices[i] = i;
554  }
555  detectordata->rowsWithConstrictions = new vector<int>();
556  detectordata->blockedAfterrow = new vector<int>();
557 
558  /* create hash tables */
559  SCIP_CALL( indexmapCreate(scip, &detectordata->indexmap, nconss, nvars) );
560 
561  return SCIP_OKAY;
562 }
563 
565 static
566 SCIP_RETCODE initData(
567  SCIP* scip,
568  gcg::Seeed* seeed,
569  DEC_DETECTORDATA* detectordata
570  )
571 {
572  int i;
573  int nvars;
574  int nconss;
575 
576  assert(seeed != NULL);
577  assert(detectordata != NULL);
578 
579  nvars = seeed->getNOpenvars();
580  nconss = seeed->getNOpenconss();
581  detectordata->maxblocks = MIN(nconss, detectordata->maxblocks);
582 
583  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->ibegin, nconss) );
584  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->iend, nconss) );
585  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jbegin, nvars) );
586  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jend, nvars) );
587  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jmin, nconss) );
588  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->jmax, nconss) );
589  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->minV, nconss-1) );
590  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->width, nconss) );
591  SCIP_CALL( SCIPallocMemoryArray(scip, &detectordata->hashmapindices, (size_t)MAX(nvars, nconss) + 1) );
592  for( i = 0; i < MAX(nvars, nconss)+1; ++i )
593  {
594  detectordata->hashmapindices[i] = i;
595  }
596  detectordata->rowsWithConstrictions = new vector<int>();
597  detectordata->blockedAfterrow = new vector<int>();
598 
599  /* create hash tables */
600  SCIP_CALL( indexmapCreate(scip, &detectordata->indexmap, nconss, nvars) );
601 
602  return SCIP_OKAY;
603 }
604 
606 static
607 SCIP_RETCODE freeData(
608  SCIP* scip,
609  DEC_DETECTORDATA* detectordata
610  )
611 {
612  assert(scip != NULL);
613  assert(detectordata != NULL);
614 
615  indexmapFree(scip, &detectordata->indexmap);
616 
617  /* delete vectors */
618  delete detectordata->rowsWithConstrictions;
619  detectordata->rowsWithConstrictions = NULL;
620  delete detectordata->blockedAfterrow;
621  detectordata->blockedAfterrow = NULL;
622 
623  if (detectordata->constoblock != NULL)
624  {
625  SCIPhashmapFree(&detectordata->constoblock);
626  detectordata->constoblock = NULL;
627  }
628 
629  SCIPfreeMemoryArray(scip, &detectordata->ibegin);
630  SCIPfreeMemoryArray(scip, &detectordata->iend);
631  SCIPfreeMemoryArray(scip, &detectordata->jbegin);
632  SCIPfreeMemoryArray(scip, &detectordata->jend);
633  SCIPfreeMemoryArray(scip, &detectordata->jmin);
634  SCIPfreeMemoryArray(scip, &detectordata->jmax);
635  SCIPfreeMemoryArray(scip, &detectordata->minV);
636  SCIPfreeMemoryArray(scip, &detectordata->width);
637  SCIPfreeMemoryArray(scip, &detectordata->hashmapindices);
638 
639  return SCIP_OKAY;
640 }
641 
658 static
659 SCIP_RETCODE createRowindexList(
660  SCIP* scip,
661  DEC_DETECTORDATA* detectordata,
662  SCIP_HASHMAP* indexcons,
663  SCIP_HASHMAP* varindex,
664  vector<vector<int> > &rowindices
665  )
666 {
667  int i;
668  int nconss = SCIPgetNConss(scip);
669 
670  for( i = 0; i < nconss; ++i )
671  {
672  int j;
673  SCIP_CONS* cons;
674  SCIP_VAR** vars;
675  int nvars;
676  int* probindices = NULL;
677  vector<int> rowindices_row;
678  int* hashmapindex = &detectordata->hashmapindices[(size_t)i+1];
679 
680  cons = (SCIP_CONS*) SCIPhashmapGetImage(indexcons, (void*) hashmapindex);
681  nvars = GCGconsGetNVars(scip, cons);
682 
683  SCIP_CALL( SCIPallocMemoryArray(scip, &vars, nvars) );
684  SCIP_CALL( GCGconsGetVars(scip, cons, vars, nvars) );
685  SCIP_CALL( SCIPallocMemoryArray(scip, &probindices, nvars) );
686 
687  /* fill the array with the indices of the variables of the current constraint */
688  for( j = 0; j < nvars; ++j )
689  {
690  probindices[j] = *(int*) SCIPhashmapGetImage(varindex, SCIPvarGetProbvar(vars[j]));
691  }
692 
693  std::sort(probindices, probindices+nvars);
694 
695  /* store a copy of the elements of probindices in the vector rowindices_row */
696  for( j = 0; j < nvars; ++j )
697  {
698  rowindices_row.push_back(probindices[j]);
699  }
700 
701  SCIPfreeMemoryArray(scip, &probindices);
702  SCIPfreeMemoryArray(scip, &vars);
703 
704  /* add rowindices_row to the vector rowindices */
705  rowindices.push_back(rowindices_row);
706  rowindices_row.clear();
707  }
708 
709  return SCIP_OKAY;
710 }
711 
728 static
729 SCIP_RETCODE createRowindexList(
730  SCIP* scip,
731  gcg::Seeed* seeed,
732  gcg::Seeedpool* seeedpool,
733  DEC_DETECTORDATA* detectordata,
734  SCIP_HASHMAP* indexcons,
735  SCIP_HASHMAP* varindex,
736  vector<vector<int> > &rowindices
737  )
738 {
739  int i;
740  int nconss = seeed->getNOpenconss();
741 
742  for( i = 0; i < nconss; ++i )
743  {
744  int j;
745  int ncurrvars;
746  int* probindices = NULL;
747  int cons;
748  vector<int> rowindices_row;
749  int* hashmapindex = &detectordata->hashmapindices[(size_t)i+1];
750 
751  cons = (int)(size_t) SCIPhashmapGetImage(indexcons, (void*) hashmapindex);
752  ncurrvars = seeedpool->getNVarsForCons(cons);
753 
754  SCIP_CALL( SCIPallocMemoryArray(scip, &probindices, ncurrvars) );
755 
756  /* fill the array with the indices of the variables of the current constraint */
757  for( j = 0; j < ncurrvars; ++j )
758  {
759  if(!seeed->isVarOpenvar(seeedpool->getVarsForCons(cons)[j]))
760  continue;
761  probindices[j] = *(int*) SCIPhashmapGetImage(varindex, (void*)(size_t)seeedpool->getVarsForCons(cons)[j]);
762  }
763 
764  std::sort(probindices, probindices+ncurrvars);
765 
766  /* store a copy of the elements of probindices in the vector rowindices_row */
767  for( j = 0; j < ncurrvars; ++j )
768  {
769  if(!seeed->isVarOpenvar(seeedpool->getVarsForCons(cons)[j]))
770  continue;
771  rowindices_row.push_back(probindices[j]);
772  }
773 
774  SCIPfreeMemoryArray(scip, &probindices);
775 
776  /* add rowindices_row to the vector rowindices */
777  rowindices.push_back(rowindices_row);
778  rowindices_row.clear();
779  }
780 
781  return SCIP_OKAY;
782 }
783 
803 static
805  SCIP* scip,
806  vector<vector<int> > &rowindices,
807  vector<vector<int> > &columnindices
808 )
809 {
810  int position;
811  int nvars;
812  int i;
813  vector<vector<int> >::iterator outer;
814  vector<int>::iterator inner;
815  nvars = SCIPgetNVars(scip);
816 
817  vector<vector<int> > columnindices_array((size_t)nvars);
818 
819  for( outer = rowindices.begin(), i = 1; outer != rowindices.end(); ++outer, ++i )
820  {
821  for( inner = outer->begin(); inner != outer->end(); ++inner )
822  {
823  position = (*inner)-1;
824  columnindices_array[(size_t) position].push_back(i);
825  }
826  }
827 
828  /* create a columnindices vector instead of an array */
829  for( i = 0; i < nvars; ++i )
830  {
831  columnindices.push_back(columnindices_array[(size_t)i]);
832  }
833 
834  return SCIP_OKAY;
835 }
836 
856 static
858  SCIP* scip,
859  gcg::Seeed* seeed,
860  vector<vector<int> > &rowindices,
861  vector<vector<int> > &columnindices
862 )
863 {
864  int position;
865  int nvars;
866  int i;
867  vector<vector<int> >::iterator outer;
868  vector<int>::iterator inner;
869  nvars = seeed->getNOpenvars();
870 
871  vector<vector<int> > columnindices_array((size_t)nvars);
872 
873  for( outer = rowindices.begin(), i = 1; outer != rowindices.end(); ++outer, ++i )
874  {
875  for( inner = outer->begin(); inner != outer->end(); ++inner )
876  {
877  position = (*inner)-1;
878  columnindices_array[(size_t) position].push_back(i);
879  }
880  }
881 
882  /* create a columnindices vector instead of an array */
883  for( i = 0; i < nvars; ++i )
884  {
885  columnindices.push_back(columnindices_array[(size_t)i]);
886  }
887 
888  return SCIP_OKAY;
889 }
890 
891 
897 static
898 vector<int> rowOrdering(
899  SCIP* scip,
900  vector<vector<int> > &columnindices,
901  int nrows
902  )
903 {
904  vector<int> roworder;
905  vector<int> new_roworder;
906  vector<vector<int> >::reverse_iterator it1;
907  vector<int>::reverse_iterator it2;
908 
909  /* create a vector for the order of the rows ( 1 2 3 ... nrows ) */
910  roworder = SCIPvectorCreateInt(1, nrows);
911  new_roworder = SCIPvectorCreateInt(1, nrows);
912 
913  /* first from back to front */
914  for( it1 = columnindices.rbegin(); it1 != columnindices.rend(); ++it1 )
915  {
916 
917  /* second from back to front */
918  for( it2 = it1->rbegin(); it2 != it1->rend(); ++it2 )
919  {
920  vector<int>::iterator tmp;
921 
922  tmp = std::find(new_roworder.begin(), new_roworder.end(), *it2); /*lint !e864*/
923  std::rotate(new_roworder.begin(), tmp, tmp+1); /*lint !e1702 !e747*/
924  }
925  roworder = new_roworder;
926  }
927 
928  return roworder;
929 }
930 
936 static
937 SCIP_RETCODE formIndexArray(
938  int* begin,
939  int* end,
940  vector<vector<int> > &indices
941  )
942 {
943  vector<vector<int> >::iterator it1;
944  int i;
945  assert(begin != NULL && end != NULL);
946  for( it1 = indices.begin(), i = 0; it1 != indices.end(); ++it1, ++i )
947  {
948  /* case: vector not empty */
949  if( !it1->empty() )
950  {
951  begin[i] = it1->front();
952  end[i] = it1->back();
953  }
954  /* case: vector empty */
955  else
956  {
957  begin[i] = 0;
958  end[i] = 0;
959  }
960  }
961  return SCIP_OKAY;
962 }
963 
964 
966 static
967 SCIP_Bool arraysAreEqual(
968  int* new_array,
969  int* old_array,
970  int num_elements
971  )
972 {
973  int i;
974  for( i = 0; i < num_elements; ++i )
975  {
976  if( new_array[i] != old_array[i] )
977  {
978  return FALSE;
979  }
980  }
981  /* case: all entries of old and new are equal */
982  return TRUE;
983 }
984 
988 static
990  SCIP* scip,
991  DEC_DETECTORDATA* detectordata,
992  INDEXMAP* inputmap,
993  INDEXMAP* outputmap
994  )
995 {
996  vector<int> roworder;
997  vector<int> columnorder;
998  vector<vector<int> > rowindices;
999  vector<vector<int> > columnindices;
1000  vector<int>::iterator it1;
1001  int nvars;
1002  int ncons;
1003  size_t i;
1004  int position;
1005  int* hashmapindex;
1006 
1007  SCIPdebugMessage("Entering rankOrderClusteringIteration\n");
1008 
1009  assert(scip != NULL);
1010  assert(detectordata != NULL);
1011  nvars = SCIPgetNVars(scip);
1012  ncons = SCIPgetNConss(scip);
1013 
1014  /* create the vectors containing the positions of nonzero entries; row and column ordering */
1015  SCIP_CALL( createRowindexList(scip, detectordata, inputmap->indexcons, inputmap->varindex, rowindices) );
1016  SCIP_CALL( createColumnindexList(scip, rowindices, columnindices) );
1017 
1018  roworder = rowOrdering(scip, columnindices, ncons);
1019  SCIPvectorRearrange(rowindices, roworder);
1020 
1021  columnorder = rowOrdering(scip, rowindices, nvars);
1022  SCIPvectorRearrange(columnindices, columnorder);
1023 
1024  /* consindex and indexcons */
1025  for( it1 = roworder.begin(), i = 0; it1 != roworder.end() && i < (size_t) ncons; ++i,++it1 )
1026  {
1027  SCIP_CONS* cons;
1028  position = *it1;
1029  hashmapindex = &detectordata->hashmapindices[position];
1030 
1031  cons = (SCIP_CONS*) SCIPhashmapGetImage(inputmap->indexcons, (void*) hashmapindex);
1032  assert ( cons != NULL);
1033 
1034  /* consindex */
1035  hashmapindex = &detectordata->hashmapindices[i+1];
1036  assert( SCIPhashmapExists(outputmap->consindex, (void*) cons));
1037  assert(*hashmapindex <= ncons);
1038  SCIP_CALL( SCIPhashmapSetImage(outputmap->consindex, (void*) cons, (void*) hashmapindex) );
1039 
1040  /* indexcons */
1041  assert( SCIPhashmapExists(outputmap->indexcons, (void*) hashmapindex ));
1042  SCIP_CALL( SCIPhashmapSetImage(outputmap->indexcons, (void*) hashmapindex, cons) );
1043  }
1044 
1045  /* varindex and indexvar */
1046  for( it1 = columnorder.begin(), i = 0; it1 != columnorder.end() &&i < (size_t) nvars; ++i, ++it1 )
1047  {
1048  SCIP_VAR* var;
1049  position = *it1;
1050  hashmapindex = &detectordata->hashmapindices[position];
1051  var = (SCIP_VAR*) SCIPhashmapGetImage(inputmap->indexvar, (void*) hashmapindex);
1052  assert ( var != NULL);
1053  assert(*hashmapindex <= nvars);
1054  /* varindex */
1055  hashmapindex = &detectordata->hashmapindices[i+1];
1056  assert( SCIPhashmapExists(outputmap->varindex, (void*) var) );
1057  SCIP_CALL( SCIPhashmapSetImage(outputmap->varindex, (void*) var, (void*) hashmapindex) );
1058 
1059  /* indexvar */
1060  assert( SCIPhashmapExists(outputmap->indexvar, (void*) hashmapindex ));
1061  SCIP_CALL( SCIPhashmapSetImage(outputmap->indexvar, (void*) hashmapindex, var) );
1062  }
1063 
1064  return SCIP_OKAY;
1065 }
1066 
1070 static
1072  SCIP* scip,
1073  gcg::Seeed* seeed,
1074  gcg::Seeedpool* seeedpool,
1075  DEC_DETECTORDATA* detectordata,
1076  INDEXMAP* inputmap,
1077  INDEXMAP* outputmap
1078  )
1079 {
1080  vector<int> roworder;
1081  vector<int> columnorder;
1082  vector<vector<int> > rowindices;
1083  vector<vector<int> > columnindices;
1084  vector<int>::iterator it1;
1085  int nvars;
1086  int ncons;
1087  size_t i;
1088  int position;
1089  int* hashmapindex;
1090 
1091  SCIPdebugMessage("Entering rankOrderClusteringIteration\n");
1092 
1093  assert(scip != NULL);
1094  assert(detectordata != NULL);
1095  nvars = seeed->getNOpenvars();
1096  ncons = seeed->getNOpenconss();
1097 
1098  /* create the vectors containing the positions of nonzero entries; row and column ordering */
1099  SCIP_CALL( createRowindexList(scip, seeed, seeedpool, detectordata, inputmap->indexcons, inputmap->varindex, rowindices) );
1100  SCIP_CALL( createColumnindexList(scip, seeed, rowindices, columnindices) );
1101 
1102  roworder = rowOrdering(scip, columnindices, ncons);
1103  SCIPvectorRearrange(rowindices, roworder);
1104 
1105  columnorder = rowOrdering(scip, rowindices, nvars);
1106  SCIPvectorRearrange(columnindices, columnorder);
1107 
1108  /* consindex and indexcons */
1109  for( it1 = roworder.begin(), i = 0; it1 != roworder.end() && i < (size_t) ncons; ++i,++it1 )
1110  {
1111  int cons;
1112  position = *it1;
1113  hashmapindex = &detectordata->hashmapindices[position];
1114 
1115  cons = (int)(size_t) SCIPhashmapGetImage(inputmap->indexcons, (void*) hashmapindex);
1116 
1117  /* consindex */
1118  hashmapindex = &detectordata->hashmapindices[i+1];
1119  assert( SCIPhashmapExists(outputmap->consindex, (void*)(size_t) cons));
1120  assert(*hashmapindex <= ncons);
1121  SCIP_CALL( SCIPhashmapSetImage(outputmap->consindex, (void*)(size_t) cons, (void*) hashmapindex) );
1122 
1123  /* indexcons */
1124  assert( SCIPhashmapExists(outputmap->indexcons, (void*) hashmapindex ));
1125  SCIP_CALL( SCIPhashmapSetImage(outputmap->indexcons, (void*) hashmapindex, (void*)(size_t) cons) );
1126  }
1127 
1128  /* varindex and indexvar */
1129  for( it1 = columnorder.begin(), i = 0; it1 != columnorder.end() &&i < (size_t) nvars; ++i, ++it1 )
1130  {
1131  int var;
1132  position = *it1;
1133  hashmapindex = &detectordata->hashmapindices[position];
1134  var = (int)(size_t) SCIPhashmapGetImage(inputmap->indexvar, (void*) hashmapindex);
1135  assert(*hashmapindex <= nvars);
1136  /* varindex */
1137  hashmapindex = &detectordata->hashmapindices[i+1];
1138  assert( SCIPhashmapExists(outputmap->varindex, (void*)(size_t) var) );
1139  SCIP_CALL( SCIPhashmapSetImage(outputmap->varindex, (void*)(size_t) var, (void*) hashmapindex) );
1140 
1141  /* indexvar */
1142  assert( SCIPhashmapExists(outputmap->indexvar, (void*) hashmapindex ));
1143  SCIP_CALL( SCIPhashmapSetImage(outputmap->indexvar, (void*) hashmapindex, (void*)(size_t) var) );
1144  }
1145 
1146  return SCIP_OKAY;
1147 }
1148 
1149 static
1150 SCIP_RETCODE rankOrderClustering(
1151  SCIP* scip,
1152  DEC_DETECTORDATA* detectordata,
1153  int max_iterations,
1154  int* iterations
1155  )
1156 {
1157  int i;
1158  INDEXMAP* indexmap_permuted;
1159  vector<vector<int> > rowindices;
1160  vector<vector<int> > columnindices;
1161  int* ibegin_permuted = NULL;
1162  int* iend_permuted = NULL;
1163  int* jbegin_permuted = NULL;
1164  int* jend_permuted = NULL;
1165  assert(scip != NULL);
1166  assert(detectordata != NULL);
1167 
1168  int nvars = SCIPgetNVars(scip);
1169  int nconss = SCIPgetNConss(scip);
1170 
1171  if( iterations != NULL )
1172  *iterations = -1;
1173 
1174  if( max_iterations <= 0 )
1175  {
1176  return SCIP_OKAY;
1177  }
1178 
1179 
1180 
1181  SCIP_CALL( indexmapCreate(scip, &indexmap_permuted, nconss, nvars) );
1182  SCIP_CALL( SCIPallocMemoryArray(scip, &ibegin_permuted, nconss) );
1183  SCIP_CALL( SCIPallocMemoryArray(scip, &iend_permuted, nconss) );
1184  SCIP_CALL( SCIPallocMemoryArray(scip, &jbegin_permuted, nvars) );
1185  SCIP_CALL( SCIPallocMemoryArray(scip, &jend_permuted, nvars) );
1186  SCIP_CALL( indexmapInit(indexmap_permuted, SCIPgetVars(scip), nvars, SCIPgetConss(scip), nconss, detectordata->hashmapindices) );
1187  i = 0;
1188 
1189  do
1190  {
1191  ++i;
1192  /* not more than max_iterations loops. no iteration limit for max_iterations == -1 */
1193  if(i > max_iterations && max_iterations != -1)
1194  break;
1195 
1196  SCIPdebugMessage("Iteration # %i of ROC2\n", i);
1197  SCIP_CALL( rankOrderClusteringIteration(scip, detectordata, detectordata->indexmap, indexmap_permuted) );
1198 
1199  /* form the new index arrays after the permutation */
1200  SCIP_CALL( createRowindexList(scip, detectordata, indexmap_permuted->indexcons, indexmap_permuted->varindex, rowindices) );
1201  SCIP_CALL( createColumnindexList(scip, rowindices, columnindices) );
1202  SCIP_CALL( formIndexArray(ibegin_permuted, iend_permuted, rowindices) );
1203  SCIP_CALL( formIndexArray(jbegin_permuted, jend_permuted, columnindices) );
1204  rowindices.clear();
1205  columnindices.clear();
1206 
1207  /* switch between index arrays containing new and old indices */
1208  swap( detectordata->ibegin, ibegin_permuted);
1209  swap( detectordata->iend, iend_permuted);
1210  swap( detectordata->jbegin, jbegin_permuted);
1211  swap( detectordata->jend, jend_permuted);
1212 
1213  /* switch between hash maps containing new and old indices */
1214  swap(detectordata->indexmap, indexmap_permuted);
1215  }
1216  /* while Index Arrays change */
1217  while( ! (arraysAreEqual(detectordata->ibegin, ibegin_permuted, nconss )
1218  && arraysAreEqual(detectordata->iend, iend_permuted, nconss)
1219  && arraysAreEqual(detectordata->jbegin, jbegin_permuted, nvars)
1220  && arraysAreEqual(detectordata->jend, jend_permuted, nvars)));
1221 
1222  indexmapFree(scip, &indexmap_permuted);
1223  SCIPfreeMemoryArray(scip, &ibegin_permuted);
1224  SCIPfreeMemoryArray(scip, &iend_permuted);
1225  SCIPfreeMemoryArray(scip, &jbegin_permuted);
1226  SCIPfreeMemoryArray(scip, &jend_permuted);
1227 
1228 
1229  if( iterations != NULL )
1230  *iterations = i-1;
1231 
1232  return SCIP_OKAY;
1233 }
1234 
1235 static
1236 SCIP_RETCODE rankOrderClustering(
1237  SCIP* scip,
1238  gcg::Seeed* seeed,
1239  gcg::Seeedpool* seeedpool,
1240  DEC_DETECTORDATA* detectordata,
1241  int max_iterations,
1242  int* iterations
1243  )
1244 {
1245  int i;
1246  INDEXMAP* indexmap_permuted;
1247  vector<vector<int> > rowindices;
1248  vector<vector<int> > columnindices;
1249  int* ibegin_permuted = NULL;
1250  int* iend_permuted = NULL;
1251  int* jbegin_permuted = NULL;
1252  int* jend_permuted = NULL;
1253  assert(scip != NULL);
1254  assert(detectordata != NULL);
1255 
1256  int nvars = seeed->getNOpenvars();
1257  int nconss = seeed->getNOpenconss();
1258 
1259  if( iterations != NULL )
1260  *iterations = -1;
1261 
1262  if( max_iterations <= 0 )
1263  {
1264  return SCIP_OKAY;
1265  }
1266 
1267 
1268 
1269  SCIP_CALL( indexmapCreate(scip, &indexmap_permuted, nconss, nvars) );
1270  SCIP_CALL( SCIPallocMemoryArray(scip, &ibegin_permuted, nconss) );
1271  SCIP_CALL( SCIPallocMemoryArray(scip, &iend_permuted, nconss) );
1272  SCIP_CALL( SCIPallocMemoryArray(scip, &jbegin_permuted, nvars) );
1273  SCIP_CALL( SCIPallocMemoryArray(scip, &jend_permuted, nvars) );
1274  SCIP_CALL( indexmapInit(indexmap_permuted, seeed, detectordata->hashmapindices) );
1275  i = 0;
1276 
1277  do
1278  {
1279  ++i;
1280  /* not more than max_iterations loops. no iteration limit for max_iterations == -1 */
1281  if(i > max_iterations && max_iterations != -1)
1282  break;
1283 
1284  SCIPdebugMessage("Iteration # %i of ROC2\n", i);
1285  SCIP_CALL( rankOrderClusteringIteration(scip, seeed, seeedpool, detectordata, detectordata->indexmap, indexmap_permuted) );
1286 
1287  /* form the new index arrays after the permutation */
1288  SCIP_CALL( createRowindexList(scip, seeed, seeedpool, detectordata, indexmap_permuted->indexcons, indexmap_permuted->varindex, rowindices) );
1289  SCIP_CALL( createColumnindexList(scip, seeed, rowindices, columnindices) );
1290  SCIP_CALL( formIndexArray(ibegin_permuted, iend_permuted, rowindices) );
1291  SCIP_CALL( formIndexArray(jbegin_permuted, jend_permuted, columnindices) );
1292  rowindices.clear();
1293  columnindices.clear();
1294 
1295  /* switch between index arrays containing new and old indices */
1296  swap( detectordata->ibegin, ibegin_permuted);
1297  swap( detectordata->iend, iend_permuted);
1298  swap( detectordata->jbegin, jbegin_permuted);
1299  swap( detectordata->jend, jend_permuted);
1300 
1301  /* switch between hash maps containing new and old indices */
1302  swap(detectordata->indexmap, indexmap_permuted);
1303  }
1304  /* while Index Arrays change */
1305  while( ! (arraysAreEqual(detectordata->ibegin, ibegin_permuted, nconss )
1306  && arraysAreEqual(detectordata->iend, iend_permuted, nconss)
1307  && arraysAreEqual(detectordata->jbegin, jbegin_permuted, nvars)
1308  && arraysAreEqual(detectordata->jend, jend_permuted, nvars)));
1309 
1310  indexmapFree(scip, &indexmap_permuted);
1311  SCIPfreeMemoryArray(scip, &ibegin_permuted);
1312  SCIPfreeMemoryArray(scip, &iend_permuted);
1313  SCIPfreeMemoryArray(scip, &jbegin_permuted);
1314  SCIPfreeMemoryArray(scip, &jend_permuted);
1315 
1316 
1317  if( iterations != NULL )
1318  *iterations = i-1;
1319 
1320  return SCIP_OKAY;
1321 }
1322 
1324 static
1326  SCIP* scip,
1327  DEC_DETECTORDATA* detectordata
1328  )
1329 {
1330  /* if blocking is performed after row i+1; local minima */
1331  size_t i;
1332  size_t nconss = (size_t) SCIPgetNConss(scip);
1333  for( i = 1; i < nconss - 2; ++i )
1334  {
1335  /* is minV[i] a local minimum? < or <= ? What does make more sense? */
1336  if( detectordata->minV[i] < detectordata->minV[i-1] && detectordata->minV[i] < detectordata->minV[i+1] )
1337  {
1338  detectordata->rowsWithConstrictions->push_back((int)(i+1));
1339  }
1340  }
1341  return SCIP_OKAY;
1342 }
1343 
1345 static
1347  SCIP* scip,
1348  gcg::Seeed* seeed,
1349  DEC_DETECTORDATA* detectordata
1350  )
1351 {
1352  /* if blocking is performed after row i+1; local minima */
1353  size_t i;
1354  size_t nconss = (size_t) seeed->getNOpenconss();
1355  for( i = 1; i < nconss - 2; ++i )
1356  {
1357  /* is minV[i] a local minimum? < or <= ? What does make more sense? */
1358  if( detectordata->minV[i] < detectordata->minV[i-1] && detectordata->minV[i] < detectordata->minV[i+1] )
1359  {
1360  detectordata->rowsWithConstrictions->push_back((int)(i+1));
1361  }
1362  }
1363  return SCIP_OKAY;
1364 }
1365 
1370 static
1371 SCIP_RETCODE assignConsToBlock(
1372  SCIP* scip,
1373  DEC_DETECTORDATA* detectordata,
1374  int block,
1375  int first_cons,
1376  int last_cons
1377  )
1378 {
1379  int i;
1380 
1381  /* assign the constraints to the current block */
1382  for( i = first_cons; i <= last_cons; ++i )
1383  {
1384  int* hashmapindex = &detectordata->hashmapindices[i];
1385  SCIP_CONS* cons = (SCIP_CONS*) SCIPhashmapGetImage(detectordata->indexmap->indexcons, (void*) hashmapindex);
1386  assert(cons != NULL);
1387  /* insert cons into hash map vartoblock */
1388  assert(!SCIPhashmapExists(detectordata->constoblock, cons));
1389  SCIP_CALL( SCIPhashmapInsert(detectordata->constoblock, cons, (void*) (size_t) (detectordata->hashmapindices[block])) );
1390  }
1391  detectordata->blockedAfterrow->push_back(detectordata->hashmapindices[last_cons]);
1392  return SCIP_OKAY;
1393 }
1394 
1396 static
1397 SCIP_RETCODE assignConsToBlock(
1398  SCIP* scip,
1399  gcg::Seeed* seeed,
1400  DEC_DETECTORDATA* detectordata,
1401  int block,
1402  int first_cons,
1403  int last_cons
1404  )
1405 {
1406  int i;
1407 
1408  /* assign the constraints to the current block */
1409  for( i = first_cons; i <= last_cons; ++i )
1410  {
1411  int* hashmapindex = &detectordata->hashmapindices[i];
1412  int cons = (int)(size_t) SCIPhashmapGetImage(detectordata->indexmap->indexcons, (void*) hashmapindex);
1413  /* insert cons into hash map vartoblock */
1414  assert(!SCIPhashmapExists(detectordata->constoblock, (void*)(size_t) cons));
1415  SCIP_CALL( SCIPhashmapInsert(detectordata->constoblock, (void*)(size_t) cons, (void*) (size_t) (detectordata->hashmapindices[block])) );
1416  }
1417  detectordata->blockedAfterrow->push_back(detectordata->hashmapindices[last_cons]);
1418  return SCIP_OKAY;
1419 }
1420 
1422 static
1424  DEC_DETECTORDATA* detectordata,
1425  int from_row,
1426  int to_row
1427  )
1428 {
1429  /* some pointer arithmetic */
1430  return std::max_element(detectordata->iend + (from_row), detectordata->iend+((size_t)to_row + 1))-detectordata->iend; /*lint !e712*/
1431 }
1432 
1434 static
1436  DEC_DETECTORDATA* detectordata,
1437  int row
1438  )
1439 {
1440  return detectordata->ibegin[row-1];
1441 }
1442 
1447 static
1449  DEC_DETECTORDATA* detectordata,
1450  int prev_block_first_row,
1451  int prev_block_last_row,
1452  int block_at_row
1453  )
1454 {
1455  int last_column_prev_block;
1456  int first_column_current_block;
1457 
1458  /* if the function is called for the first block, the blocking is always valid */
1459  if( prev_block_last_row == 0 )
1460  {
1461  return TRUE;
1462  }
1463  last_column_prev_block = getMaxColIndex(detectordata, prev_block_first_row, prev_block_last_row);
1464  first_column_current_block = getMinColIndex(detectordata, block_at_row);
1465  return ( first_column_current_block > last_column_prev_block ? TRUE : FALSE);
1466 }
1467 
1472 static
1473 vector<int>::iterator findBlockingCandidate(
1474  vector<int>::iterator it_constrictions,
1475  vector<int>* it_vector,
1476  int min_block_size,
1477  int prev_block_last_row
1478  )
1479 {
1480  for( ;; )
1481  {
1482  /* end of the vector? */
1483  if( it_constrictions == it_vector->end() )
1484  {
1485  return it_constrictions;
1486  }
1487  /* does a blocking to the next row forming a constriction comprise more rows than min_block_size */
1488  if( (*it_constrictions - prev_block_last_row) >= min_block_size )
1489  {
1490  return it_constrictions;
1491  }
1492  /* advance iterator to next element */
1493  ++it_constrictions;
1494  }
1495 }
1496 
1501 static
1502 vector<int>::iterator nextRowToBlockAt(
1503  DEC_DETECTORDATA* detectordata,
1504  vector<int>::iterator it_constrictions,
1505  vector<int>* it_vector,
1506  int min_block_size,
1507  int prev_block_first_row,
1508  int prev_block_last_row
1509  )
1510 {
1511 
1512  /* end of the constriction vector? */
1513  if( it_constrictions == it_vector->end() )
1514  {
1515  return it_constrictions;
1516  }
1517 
1518  for( ;; )
1519  {
1520  /* find a blocking candidate */
1521  it_constrictions = findBlockingCandidate(it_constrictions, it_vector, min_block_size, prev_block_last_row);
1522  /* case: no candidate found */
1523  if( it_constrictions == it_vector->end() )
1524  {
1525  break;
1526  }
1527  /* case: candidate found */
1528  else
1529  {
1530  /* valid blocking */
1531  if( isValidBlocking(detectordata, prev_block_first_row, prev_block_last_row, *it_constrictions) )
1532  {
1533  break;
1534  }
1535  /* invalid blocking */
1536  else
1537  {
1538  ++it_constrictions;
1539  }
1540  }
1541  }
1542  return it_constrictions;
1543 }
1544 
1546 static
1548  DEC_DETECTORDATA* detectordata
1549  )
1550 {
1551  int nblockingtypes;
1552  int nblockingspertype;
1553 
1554  nblockingtypes = 0;
1555  /* get the number of enabled blocking types */
1556  if( detectordata->dynamicblocking )
1557  {
1558  ++nblockingtypes;
1559  }
1560  if( detectordata->staticblocking )
1561  {
1562  ++nblockingtypes;
1563  }
1564  if( detectordata->blockingassoonaspossible )
1565  {
1566  ++nblockingtypes;
1567  }
1568 
1569  /* get the number of blockings per blocking type */
1570  if( detectordata->multipledecomps )
1571  {
1572  nblockingspertype = detectordata->maxblocks - detectordata->minblocks + 1;
1573  }
1574  else
1575  {
1576  nblockingspertype = 1;
1577  }
1578 
1579  return nblockingtypes * nblockingspertype;
1580 }
1581 
1583 static
1585  DEC_DETECTORDATA* detectordata,
1586  SCIP_RESULT* result
1587  )
1588 {
1589  /* maxblocks < nRelevantsCons? */
1590 
1591  /* desired blocks <= maxblocks? */
1592 
1593  /* is minblocks <= maxblocks? */
1594  if( detectordata->multipledecomps )
1595  {
1596  if( detectordata->minblocks > detectordata->maxblocks )
1597  {
1598  SCIPerrorMessage("minblocks > maxblocks. Setting minblocks = maxblocks.\n");
1599  detectordata->minblocks = detectordata->maxblocks;
1600  }
1601  }
1602 
1603  /* is at least one blocking type enabled? */
1604  if( ! detectordata->blockingassoonaspossible && ! detectordata->staticblocking &&! detectordata->dynamicblocking )
1605  {
1606  SCIPerrorMessage("No blocking type enabled, cannot perform blocking.\n");
1607  *result = SCIP_DIDNOTRUN;
1608  }
1609 }
1610 
1612 static
1613 SCIP_RETCODE blockingDynamic(
1614  SCIP* scip,
1615  DEC_DETECTORDATA* detectordata,
1616  int tau,
1617  int nvars
1618  )
1619 { /*lint -e715*/
1620  int block;
1621  int prev_block_first_row;
1622  int prev_block_last_row;
1623  int min_block_size;
1624  /* notation: i=current block; im1=i-1=previous block; ip1=i+1=next block */
1625 #ifdef SCIP_DEBUG
1626  int max_col_index_im1 = 0;
1627 #endif
1628  vector<int>::iterator it1;
1629  /* debug */
1630  SCIPdebugMessage("Starting Blocking...\n");
1631  SCIPdebugMessage("Max blocks: %i\n", detectordata->maxblocks);
1632  block = 1;
1633  prev_block_first_row = 0;
1634  prev_block_last_row = 0;
1635 
1636  assert(tau > 0);
1637  min_block_size = (int) SCIPround(scip, ((SCIP_Real)SCIPgetNConss(scip)) / (2.0 * tau ));
1638  it1 = detectordata->rowsWithConstrictions->begin();
1639 
1640  for( it1 = nextRowToBlockAt(detectordata, it1, detectordata->rowsWithConstrictions, min_block_size, prev_block_first_row, prev_block_last_row);
1641  it1 != detectordata->rowsWithConstrictions->end() && block < detectordata->maxblocks;
1642  it1 = nextRowToBlockAt(detectordata, it1, detectordata->rowsWithConstrictions, min_block_size, prev_block_first_row, prev_block_last_row) )
1643  {
1644  int current_row = * it1;
1645 #ifdef SCIP_DEBUG
1646  int max_col_index_i = getMaxColIndex(detectordata, prev_block_last_row + 1, current_row);
1647  int min_col_index_ip1 = getMinColIndex(detectordata, current_row + 1);
1648  SCIPdebugMessage("vars in block: %i - %i, linking vars: %i - %i\n", max_col_index_im1+1, max_col_index_i, min_col_index_ip1, max_col_index_i);
1649 #endif
1650  /* assign the variables and constraints to block */
1651  SCIP_CALL( assignConsToBlock(scip, detectordata, block, prev_block_last_row + 1, current_row) );
1652  /* update variables in the while loop */
1653  prev_block_first_row = prev_block_last_row + 1;
1654  prev_block_last_row = current_row;
1655  ++block;
1656 
1657 #ifdef SCIP_DEBUG
1658  max_col_index_im1 = max_col_index_i;
1659 #endif
1660  }
1661  /* assign the remaining (< M/2tau) cons and vars to the last block; no new linking vars are added */
1662 #ifdef SCIP_DEBUG
1663  SCIPdebugMessage("last time: vars in block: %i - %i, linking vars: %i - %i\n", max_col_index_im1+1, nvars, nvars+1, nvars);
1664 #endif
1665  SCIP_CALL( assignConsToBlock(scip, detectordata, block, prev_block_last_row + 1, SCIPgetNConss(scip)) );
1666  detectordata->blockedAfterrow->pop_back();
1667 
1668  detectordata->blocks = block;
1669 
1670  /* debug plot the blocking plot for [i=1:2:1] 'test.dat' every :::i::i lt i pt 5 */
1671 #ifdef WRITEALLOUTPUT
1672  {
1673  char filename[256];
1674  char paramfile[256];
1675 
1676  sprintf(filename, "%s_dynamic_minV", getProbNameWithoutPath(scip));
1677  sprintf(paramfile, "%s_dynamic.params", getProbNameWithoutPath(scip));
1678  plotMinV(scip, detectordata, filename);
1679  }
1680 #endif
1681 
1682  return SCIP_OKAY;
1683 }
1684 
1686 static
1687 SCIP_RETCODE blockingDynamic(
1688  SCIP* scip,
1689  gcg::Seeed* seeed,
1690  DEC_DETECTORDATA* detectordata,
1691  int tau,
1692  int nvars
1693  )
1694 { /*lint -e715*/
1695  int block;
1696  int prev_block_first_row;
1697  int prev_block_last_row;
1698  int min_block_size;
1699  /* notation: i=current block; im1=i-1=previous block; ip1=i+1=next block */
1700 #ifdef SCIP_DEBUG
1701  int max_col_index_im1 = 0;
1702 #endif
1703  vector<int>::iterator it1;
1704  /* debug */
1705  SCIPdebugMessage("Starting Blocking...\n");
1706  SCIPdebugMessage("Max blocks: %i\n", detectordata->maxblocks);
1707  block = 1;
1708  prev_block_first_row = 0;
1709  prev_block_last_row = 0;
1710 
1711  assert(tau > 0);
1712  min_block_size = (int) SCIPround(scip, ((SCIP_Real)seeed->getNOpenconss()) / (2.0 * tau ));
1713  it1 = detectordata->rowsWithConstrictions->begin();
1714 
1715  for( it1 = nextRowToBlockAt(detectordata, it1, detectordata->rowsWithConstrictions, min_block_size, prev_block_first_row, prev_block_last_row);
1716  it1 != detectordata->rowsWithConstrictions->end() && block < detectordata->maxblocks;
1717  it1 = nextRowToBlockAt(detectordata, it1, detectordata->rowsWithConstrictions, min_block_size, prev_block_first_row, prev_block_last_row) )
1718  {
1719  int current_row = * it1;
1720 #ifdef SCIP_DEBUG
1721  int max_col_index_i = getMaxColIndex(detectordata, prev_block_last_row + 1, current_row);
1722  int min_col_index_ip1 = getMinColIndex(detectordata, current_row + 1);
1723  SCIPdebugMessage("vars in block: %i - %i, linking vars: %i - %i\n", max_col_index_im1+1, max_col_index_i, min_col_index_ip1, max_col_index_i);
1724 #endif
1725  /* assign the variables and constraints to block */
1726  SCIP_CALL( assignConsToBlock(scip, seeed, detectordata, block, prev_block_last_row + 1, current_row) );
1727  /* update variables in the while loop */
1728  prev_block_first_row = prev_block_last_row + 1;
1729  prev_block_last_row = current_row;
1730  ++block;
1731 
1732 #ifdef SCIP_DEBUG
1733  max_col_index_im1 = max_col_index_i;
1734 #endif
1735  }
1736  /* assign the remaining (< M/2tau) cons and vars to the last block; no new linking vars are added */
1737 #ifdef SCIP_DEBUG
1738  SCIPdebugMessage("last time: vars in block: %i - %i, linking vars: %i - %i\n", max_col_index_im1+1, nvars, nvars+1, nvars);
1739 #endif
1740  SCIP_CALL( assignConsToBlock(scip, seeed, detectordata, block, prev_block_last_row + 1, SCIPgetNConss(scip)) );
1741  detectordata->blockedAfterrow->pop_back();
1742 
1743  detectordata->blocks = block;
1744 
1745  /* debug plot the blocking plot for [i=1:2:1] 'test.dat' every :::i::i lt i pt 5 */
1746 #ifdef WRITEALLOUTPUT
1747  {
1748  char filename[256];
1749  char paramfile[256];
1750 
1751  sprintf(filename, "%s_dynamic_minV", getProbNameWithoutPath(scip));
1752  sprintf(paramfile, "%s_dynamic.params", getProbNameWithoutPath(scip));
1753  plotMinV(scip, detectordata, filename);
1754  }
1755 #endif
1756 
1757  return SCIP_OKAY;
1758 }
1759 
1761 static
1762 SCIP_RETCODE blockingStatic(
1763  SCIP* scip,
1764  DEC_DETECTORDATA* detectordata
1765  )
1766 {
1767  int nblocks;
1768  int block;
1769  int prev_block_last_row;
1770  int current_row;
1771  int nconss;
1772  /* notation: i=current block; im1=i-1=previous block; ip1=i+1=next block */
1773 
1774  assert(scip != NULL);
1775  assert(detectordata != NULL);
1776  nconss = SCIPgetNConss(scip);
1777  nblocks = nconss/detectordata->nconssperblock;
1778  prev_block_last_row = 0;
1779  current_row = 0;
1780 
1781  /* blocks 1 to (desired_blocks-1) */
1782  for( block = 1; block <= nblocks; ++block )
1783  {
1784  current_row += detectordata->nconssperblock;
1785  SCIPdebugMessage("Blocking from %d to %d in block %d",prev_block_last_row + 1, current_row, block);
1786  SCIP_CALL( assignConsToBlock(scip, detectordata, block, prev_block_last_row + 1, current_row) );
1787 
1788 
1789  prev_block_last_row = current_row;
1790  }
1791  /* last block */
1792  /* assign the remaining cons and vars to the last block; no new linking vars are added */
1793  if( prev_block_last_row + 1 <= nconss)
1794  {
1795  SCIPdebugMessage("last time: assignVarsToBlock: block, from_row, to_row: %i, %i, %i\n", block, prev_block_last_row + 1, SCIPgetNConss(scip));
1796 
1797  SCIP_CALL( assignConsToBlock(scip, detectordata, block, prev_block_last_row + 1, nconss) );
1798  detectordata->blockedAfterrow->pop_back();
1799  ++block;
1800  }
1801  detectordata->blocks = block-1;
1802 
1803 #ifdef WRITEALLOUTPUT
1804  {
1805  char filename[256];
1806  char paramfile[256];
1807 
1808  /* debug */
1809  sprintf(filename, "%s_static_minV_%i", getProbNameWithoutPath(scip), detectordata->blocks);
1810  sprintf(paramfile, "%s_static.params", getProbNameWithoutPath(scip));
1811  plotMinV(scip, detectordata, filename);
1812  }
1813 #endif
1814 
1815  return SCIP_OKAY;
1816 }
1817 
1819 static
1820 SCIP_RETCODE blockingStatic(
1821  SCIP* scip,
1822  gcg::Seeed* seeed,
1823  DEC_DETECTORDATA* detectordata
1824  )
1825 {
1826  int nblocks;
1827  int block;
1828  int prev_block_last_row;
1829  int current_row;
1830  int nconss;
1831  /* notation: i=current block; im1=i-1=previous block; ip1=i+1=next block */
1832 
1833  assert(scip != NULL);
1834  assert(detectordata != NULL);
1835  nconss = seeed->getNOpenconss();
1836  nblocks = nconss/detectordata->nconssperblock;
1837  prev_block_last_row = 0;
1838  current_row = 0;
1839 
1840  /* blocks 1 to (desired_blocks-1) */
1841  for( block = 1; block <= nblocks; ++block )
1842  {
1843  current_row += detectordata->nconssperblock;
1844  SCIPdebugMessage("Blocking from %d to %d in block %d",prev_block_last_row + 1, current_row, block);
1845  SCIP_CALL( assignConsToBlock(scip, seeed, detectordata, block, prev_block_last_row + 1, current_row) );
1846 
1847 
1848  prev_block_last_row = current_row;
1849  }
1850  /* last block */
1851  /* assign the remaining cons and vars to the last block; no new linking vars are added */
1852  if( prev_block_last_row + 1 <= nconss)
1853  {
1854  SCIPdebugMessage("last time: assignVarsToBlock: block, from_row, to_row: %i, %i, %i\n", block, prev_block_last_row + 1, SCIPgetNConss(scip));
1855 
1856  SCIP_CALL( assignConsToBlock(scip, seeed, detectordata, block, prev_block_last_row + 1, nconss) );
1857  detectordata->blockedAfterrow->pop_back();
1858  ++block;
1859  }
1860  detectordata->blocks = block-1;
1861 
1862 #ifdef WRITEALLOUTPUT
1863  {
1864  char filename[256];
1865  char paramfile[256];
1866 
1867  /* debug */
1868  sprintf(filename, "%s_static_minV_%i", getProbNameWithoutPath(scip), detectordata->blocks);
1869  sprintf(paramfile, "%s_static.params", getProbNameWithoutPath(scip));
1870  plotMinV(scip, detectordata, filename);
1871  }
1872 #endif
1873 
1874  return SCIP_OKAY;
1875 }
1876 
1877 static
1879  SCIP* scip,
1880  DEC_DETECTORDATA* detectordata,
1881  int desired_blocks,
1882  int nvars
1883 )
1884 {
1885  int block;
1886  block = 0;
1887  detectordata->blocks = block;
1888 
1889  return SCIP_OKAY;
1890 }
1891 
1893 static
1894 SCIP_RETCODE resetDetectordata(
1895  SCIP* scip,
1896  DEC_DETECTORDATA* detectordata
1897  )
1898 {
1899  if(detectordata->constoblock != NULL)
1900  {
1901  SCIP_CALL( SCIPhashmapRemoveAll(detectordata->constoblock) );
1902  }
1903  else
1904  {
1905  SCIP_CALL( SCIPhashmapCreate(&(detectordata->constoblock), SCIPblkmem(scip), SCIPgetNConss(scip)) );
1906  }
1907  return SCIP_OKAY;
1908 }
1909 
1910 static
1911 SCIP_RETCODE blocking(
1912  SCIP* scip,
1913  DEC_DETECTORDATA* detectordata,
1914  DEC_DECOMP*** decdecomps,
1915  int* ndecdecomps,
1916  int nvars,
1917  int ncons,
1918  SCIP_RESULT* result
1919  )
1920 {
1921  int tau = 1; /* desired number of blocks */
1922 
1923  assert(*ndecdecomps == 0);
1924 
1925  SCIPdebugMessage("Entering Blocking\n");
1926 
1927  /* if multiple decompositions disabled */
1928  if( detectordata->multipledecomps == FALSE )
1929  {
1930  /* if desiredblocks == 0 let the algorithm determine the desired number of blocks */
1931  if( detectordata->desiredblocks == 0 )
1932  {
1933  int n = *std::max_element(detectordata->width, detectordata->width+ncons);
1934  int v = *std::min_element(detectordata->width, detectordata->width+ncons);
1935  tau = (int)SCIPround(scip, ((SCIP_Real)nvars - v)/((SCIP_Real)n - v));
1936  SCIPdebugMessage("<n><v><tau>: <%i><%i><%i>\n", n, v, tau);
1937  if( tau > detectordata->maxblocks )
1938  {
1939  tau = detectordata->maxblocks;
1940  }
1941 
1942  SCIPdebugMessage("detectordata->enablemultipledecomps = 0. detectordata->desiredblocks == 0. Calculating tau = %i\n", tau);
1943  /* continue only if tau >= 2 */
1944  if( tau < 2 )
1945  {
1946  *result = SCIP_DIDNOTFIND;
1947  return SCIP_OKAY;
1948  }
1949  }
1950  else
1951  {
1952  tau = detectordata->desiredblocks;
1953  }
1954  }
1955 
1956  /* variant 2 */
1957  /* dynamic blocking */
1958  if( detectordata->dynamicblocking )
1959  {
1960  SCIPdebugMessage("detectordata->enableblockingdynamic = 1.\n");
1961  SCIP_CALL( rowsWithConstriction(scip, detectordata) );
1962 
1963  SCIPdebugMessage("detectordata->enablemultipledecomps = %ud.\n", detectordata->multipledecomps);
1964 
1965  if( detectordata->multipledecomps )
1966  {
1967  for( tau = detectordata->minblocks; tau <= detectordata->maxblocks; ++tau )
1968  {
1969  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *ndecdecomps);
1970  SCIP_CALL( resetDetectordata(scip, detectordata) );
1971 
1972  SCIP_CALL( blockingDynamic(scip, detectordata, tau, nvars) );
1973  if( detectordata->blocks <= 1 )
1974  continue;
1975 
1976  SCIP_CALL( DECdecompCreate(scip, &((*decdecomps)[*ndecdecomps])) );
1977  SCIP_CALL( DECfilloutDecompFromConstoblock(scip, (*decdecomps)[*ndecdecomps], detectordata->constoblock, detectordata->blocks, TRUE) );
1978  // detectordata->constoblock = NULL;
1979 
1980  (*ndecdecomps) += 1;
1981  }
1982  }
1983  else
1984  {
1985  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *ndecdecomps);
1986  SCIP_CALL( resetDetectordata(scip, detectordata) );
1987 
1988  SCIP_CALL( blockingDynamic(scip, detectordata, tau, nvars) );
1989  if( detectordata->blocks > 1 )
1990  {
1991  SCIP_CALL( DECdecompCreate(scip, &((*decdecomps)[*ndecdecomps])) );
1992  SCIP_CALL( DECfilloutDecompFromConstoblock(scip, (*decdecomps)[*ndecdecomps], detectordata->constoblock, detectordata->blocks, TRUE) );
1993  //detectordata->constoblock = NULL;
1994 
1995  (*ndecdecomps) += 1;
1996  }
1997  }
1998  }
1999 
2000  /* static blocking */
2001  SCIPdebugMessage("detectordata->staticblocking = %ud. \n", detectordata->staticblocking);
2002 
2003  if( detectordata->staticblocking )
2004  {
2005  SCIPdebugMessage("nconssperblock = %i, dec = %i\n", detectordata->nconssperblock, *ndecdecomps);
2006 
2007  SCIP_CALL( resetDetectordata(scip, detectordata) );
2008  SCIP_CALL( blockingStatic(scip, detectordata) );
2009 
2010  if( detectordata->blocks > 1 )
2011  {
2012  SCIP_CALL( DECdecompCreate(scip, &((*decdecomps)[*ndecdecomps])) );
2013  SCIP_CALL( DECfilloutDecompFromConstoblock(scip, (*decdecomps)[*ndecdecomps], detectordata->constoblock, detectordata->blocks, TRUE) );
2014  //detectordata->constoblock = NULL;
2015 
2016  (*ndecdecomps) += 1;
2017  }
2018  }
2019 
2020  /* blocking ASAP */
2021  SCIPdebugMessage("detectordata->blockingassoonaspossible = %ud. \n", detectordata->blockingassoonaspossible);
2022 
2023  if( detectordata->blockingassoonaspossible )
2024  {
2025  if( detectordata->multipledecomps )
2026  {
2027  for( tau = detectordata->minblocks; tau <= detectordata->maxblocks; ++tau )
2028  {
2029  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *ndecdecomps);
2030  SCIP_CALL( resetDetectordata(scip, detectordata) );
2031 
2032  SCIP_CALL( blockingAsSoonAsPossible(scip, detectordata, tau, nvars) );
2033  if( detectordata->blocks <= 1)
2034  continue;
2035 
2036  SCIP_CALL( DECdecompCreate(scip, &((*decdecomps)[*ndecdecomps])) );
2037  SCIP_CALL( DECfilloutDecompFromConstoblock(scip, (*decdecomps)[*ndecdecomps], detectordata->constoblock, detectordata->blocks, TRUE) );
2038  detectordata->constoblock = NULL;
2039 
2040  *ndecdecomps += 1;
2041  }
2042  }
2043  else
2044  {
2045  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *ndecdecomps);
2046  SCIP_CALL( resetDetectordata(scip, detectordata) );
2047  SCIP_CALL( blockingAsSoonAsPossible(scip, detectordata, tau, nvars) );
2048  if( detectordata->blocks > 1)
2049  {
2050  SCIP_CALL( DECdecompCreate(scip, &((*decdecomps)[*ndecdecomps])) );
2051  SCIP_CALL( DECfilloutDecompFromConstoblock(scip, (*decdecomps)[*ndecdecomps], detectordata->constoblock, detectordata->blocks, TRUE) );
2052  detectordata->constoblock = NULL;
2053 
2054  *ndecdecomps += 1;
2055  }
2056  }
2057  }
2058  return SCIP_OKAY;
2059 }
2060 
2061 static
2062 SCIP_RETCODE blocking(
2063  SCIP* scip,
2064  DEC_DETECTORDATA* detectordata,
2065  gcg::Seeed* seeed,
2066  gcg::Seeedpool* seeedpool,
2067  gcg::Seeed*** newSeeeds,
2068  int* nNewSeeeds,
2069  SCIP_Real time,
2070  SCIP_RESULT* result
2071  )
2072 {
2073  SCIP_CLOCK* temporaryClock;
2074  SCIP_CALL_ABORT( SCIPcreateClock(scip, &temporaryClock) );
2075  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2076 
2077  int tau = 1; /* desired number of blocks */
2078  int ncons = seeed->getNOpenconss();
2079  int nvars = seeed->getNOpenvars();
2080 
2081  assert(*nNewSeeeds == 0);
2082 
2083  SCIPdebugMessage("Entering Blocking\n");
2084 
2085  /* if multiple decompositions disabled */
2086  if( detectordata->multipledecomps == FALSE )
2087  {
2088  /* if desiredblocks == 0 let the algorithm determine the desired number of blocks */
2089  if( detectordata->desiredblocks == 0 )
2090  {
2091  int n = *std::max_element(detectordata->width, detectordata->width+ncons);
2092  int v = *std::min_element(detectordata->width, detectordata->width+ncons);
2093  tau = (int)SCIPround(scip, ((SCIP_Real)nvars - v)/((SCIP_Real)n - v));
2094  SCIPdebugMessage("<n><v><tau>: <%i><%i><%i>\n", n, v, tau);
2095  if( tau > detectordata->maxblocks )
2096  {
2097  tau = detectordata->maxblocks;
2098  }
2099 
2100  SCIPdebugMessage("detectordata->enablemultipledecomps = 0. detectordata->desiredblocks == 0. Calculating tau = %i\n", tau);
2101  /* continue only if tau >= 2 */
2102  if( tau < 2 )
2103  {
2104  *result = SCIP_DIDNOTFIND;
2105  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2106  SCIP_CALL_ABORT(SCIPfreeClock(scip, &temporaryClock) );
2107  return SCIP_OKAY;
2108  }
2109  }
2110  else
2111  {
2112  tau = detectordata->desiredblocks;
2113  }
2114  }
2115 
2116  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2117  SCIP_Real tempTime = SCIPclockGetTime(temporaryClock);
2118  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2119  /* variant 2 */
2120  /* dynamic blocking */
2121  if( detectordata->dynamicblocking )
2122  {
2123  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2124  SCIPdebugMessage("detectordata->enableblockingdynamic = 1.\n");
2125  SCIP_CALL( rowsWithConstriction(scip, seeed, detectordata) );
2126 
2127  SCIPdebugMessage("detectordata->enablemultipledecomps = %ud.\n", detectordata->multipledecomps);
2128 
2129  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2130  SCIP_Real tempTimeDynamicBlocking = SCIPclockGetTime(temporaryClock);
2131  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2132 
2133  if( detectordata->multipledecomps )
2134  {
2135  for( tau = detectordata->minblocks; tau <= detectordata->maxblocks; ++tau )
2136  {
2137  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2138  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *nNewSeeeds);
2139  SCIP_CALL( resetDetectordata(scip, detectordata) );
2140 
2141  SCIP_CALL( blockingDynamic(scip, seeed, detectordata, tau, nvars) );
2142  if( detectordata->blocks <= 1 )
2143  {
2144  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2145  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2146  continue;
2147  }
2148 
2149  (*newSeeeds[*nNewSeeeds]) = new gcg::Seeed(seeed);
2150  SCIP_CALL((*newSeeeds[*nNewSeeeds])->assignSeeedFromConstoblock(detectordata->constoblock, detectordata->blocks) );
2151  (*newSeeeds[*nNewSeeeds])->assignCurrentStairlinking();
2152 
2153  // detectordata->constoblock = NULL;
2154 
2155  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2156  (*newSeeeds[*nNewSeeeds])->addClockTime( time + tempTime + tempTimeDynamicBlocking + SCIPclockGetTime(temporaryClock) );
2157  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2158 
2159  (*nNewSeeeds) += 1;
2160  }
2161  }
2162  else
2163  {
2164  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2165  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *nNewSeeeds);
2166  SCIP_CALL( resetDetectordata(scip, detectordata) );
2167 
2168  SCIP_CALL( blockingDynamic(scip, seeed, detectordata, tau, nvars) );
2169  if( detectordata->blocks > 1 )
2170  {
2171  (*newSeeeds[*nNewSeeeds]) = new gcg::Seeed(seeed);
2172  SCIP_CALL((*newSeeeds[*nNewSeeeds])->assignSeeedFromConstoblock(detectordata->constoblock, detectordata->blocks) );
2173  (*newSeeeds[*nNewSeeeds])->assignCurrentStairlinking();
2174  // detectordata->constoblock = NULL;
2175 
2176  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2177  (*newSeeeds[*nNewSeeeds])->addClockTime( time + tempTime + tempTimeDynamicBlocking + SCIPclockGetTime(temporaryClock) );
2178  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2179 
2180  (*nNewSeeeds) += 1;
2181  }
2182  else
2183  {
2184  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2185  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2186  }
2187  }
2188  }
2189 
2190  /* static blocking */
2191  SCIPdebugMessage("detectordata->staticblocking = %ud. \n", detectordata->staticblocking);
2192 
2193  if( detectordata->staticblocking )
2194  {
2195  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2196  SCIPdebugMessage("nconssperblock = %i, dec = %i\n", detectordata->nconssperblock, *nNewSeeeds);
2197 
2198  SCIP_CALL( resetDetectordata(scip, detectordata) );
2199  SCIP_CALL( blockingStatic(scip, seeed, detectordata) );
2200 
2201  if( detectordata->blocks > 1 )
2202  {
2203  (*newSeeeds[*nNewSeeeds]) = new gcg::Seeed(seeed);
2204  SCIP_CALL((*newSeeeds[*nNewSeeeds])->assignSeeedFromConstoblock(detectordata->constoblock, detectordata->blocks) );
2205  (*newSeeeds[*nNewSeeeds])->assignCurrentStairlinking();
2206 // detectordata->constoblock = NULL;
2207 
2208  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2209  (*newSeeeds[*nNewSeeeds])->addClockTime( time + tempTime + SCIPclockGetTime(temporaryClock) );
2210  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2211 
2212  (*nNewSeeeds) += 1;
2213  }
2214  else
2215  {
2216  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2217  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2218  }
2219  }
2220 
2221  /* blocking ASAP */
2222  SCIPdebugMessage("detectordata->blockingassoonaspossible = %ud. \n", detectordata->blockingassoonaspossible);
2223 
2224  if( detectordata->blockingassoonaspossible )
2225  {
2226  if( detectordata->multipledecomps )
2227  {
2228  for( tau = detectordata->minblocks; tau <= detectordata->maxblocks; ++tau )
2229  {
2230  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2231  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *nNewSeeeds);
2232  SCIP_CALL( resetDetectordata(scip, detectordata) );
2233 
2234  SCIP_CALL( blockingAsSoonAsPossible(scip, detectordata, tau, nvars) );
2235  if( detectordata->blocks <= 1)
2236  {
2237  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2238  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2239  continue;
2240  }
2241 
2242  (*newSeeeds[*nNewSeeeds]) = new gcg::Seeed(seeed);
2243  SCIP_CALL((*newSeeeds[*nNewSeeeds])->assignSeeedFromConstoblock(detectordata->constoblock, detectordata->blocks) );
2244  (*newSeeeds[*nNewSeeeds])->assignCurrentStairlinking();
2245  // detectordata->constoblock = NULL;
2246 
2247  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2248  (*newSeeeds[*nNewSeeeds])->addClockTime( time + tempTime + SCIPclockGetTime(temporaryClock) );
2249  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2250 
2251  *nNewSeeeds += 1;
2252  }
2253  }
2254  else
2255  {
2256  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2257  SCIPdebugMessage("tau = %i, dec = %i\n", tau, *nNewSeeeds);
2258  SCIP_CALL( resetDetectordata(scip, detectordata) );
2259  SCIP_CALL( blockingAsSoonAsPossible(scip, detectordata, tau, nvars) );
2260  if( detectordata->blocks > 1)
2261  {
2262  (*newSeeeds[*nNewSeeeds]) = new gcg::Seeed(seeed);
2263  SCIP_CALL((*newSeeeds[*nNewSeeeds])->assignSeeedFromConstoblock(detectordata->constoblock, detectordata->blocks) );
2264  (*newSeeeds[*nNewSeeeds])->assignCurrentStairlinking();
2265  if( detectordata->constoblock != NULL )
2266  SCIPhashmapFree( &detectordata->constoblock );
2267  detectordata->constoblock = NULL;
2268 
2269  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2270  (*newSeeeds[*nNewSeeeds])->addClockTime( time + tempTime + SCIPclockGetTime(temporaryClock) );
2271  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2272 
2273  *nNewSeeeds += 1;
2274  }
2275  else
2276  {
2277  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock) );
2278  SCIP_CALL_ABORT( SCIPresetClock(scip, temporaryClock) );
2279  }
2280  }
2281  }
2282  SCIP_CALL_ABORT(SCIPfreeClock(scip, &temporaryClock) );
2283  return SCIP_OKAY;
2284 }
2285 
2286 
2288 static
2289 DEC_DECL_FREEDETECTOR(detectorFreeStairheur)
2290 {
2291  DEC_DETECTORDATA* detectordata;
2292 
2293  assert(scip != NULL);
2294 
2295  detectordata = DECdetectorGetData(detector);
2296  assert(detectordata != NULL);
2297 
2298  if (detectordata->constoblock != NULL)
2299  {
2300  SCIPhashmapFree(&detectordata->constoblock);
2301  detectordata->constoblock = NULL;
2302  }
2303 
2304  assert(strcmp(DECdetectorGetName(detector), DEC_DETECTORNAME) == 0);
2305 
2306  SCIPfreeMemory(scip, &detectordata);
2307  return SCIP_OKAY;
2308 }
2309 
2311 static
2312 DEC_DECL_INITDETECTOR(detectorInitStairheur)
2313 {
2314  DEC_DETECTORDATA* detectordata;
2315 
2316  assert(scip != NULL);
2317 
2318  detectordata = DECdetectorGetData(detector);
2319  assert(detectordata != NULL);
2320 
2321  detectordata->constoblock = NULL;
2322  detectordata->ibegin = NULL;
2323  detectordata->iend = NULL;
2324  detectordata->jbegin = NULL;
2325  detectordata->jend = NULL;
2326  detectordata->jmin = NULL;
2327  detectordata->jmax = NULL;
2328  detectordata->minV = NULL;
2329  detectordata->width = NULL;
2330  detectordata->hashmapindices = NULL;
2331  detectordata->indexmap = NULL;
2332  detectordata->rowsWithConstrictions = NULL;
2333  detectordata->blockedAfterrow = NULL;
2334 
2335  return SCIP_OKAY;
2336 }
2337 
2339 static
2340 DEC_DECL_DETECTSTRUCTURE(detectorDetectStairheur)
2341 {
2342  int i;
2343  int nconss; /* number of constraints in the problem */
2344  int nvars; /* number of variables in the problem */
2345  int ndecs;
2346  SCIP_VAR** vars;
2347  SCIP_CONS** conss;
2348  vector<vector<int> > rowindices;
2349  vector<vector<int> > columnindices;
2350 #ifdef WRITEALLOUTPUT
2351  int ROC_iterations;
2352 #endif
2353 
2354  assert(scip != NULL);
2355  assert(detectordata != NULL);
2356  assert(decdecomps != NULL);
2357  assert(ndecdecomps != NULL);
2358 
2359  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, "Detecting stairheur structure:");
2360 
2361  SCIP_CALL( initData(scip, detectordata) );
2362 
2363  checkParameterConsistency(detectordata, result);
2364  ndecs = calculateNdecompositions(detectordata);
2365  SCIPdebugMessage("%i decompositions will be created\n", ndecs);
2366  *ndecdecomps = 0;
2367 
2368  /* allocate space for output data */
2369  SCIP_CALL( SCIPallocMemoryArray(scip, decdecomps, ndecs) );
2370 
2371  nvars = SCIPgetNVars(scip);
2372  vars = SCIPgetVars(scip);
2373  nconss = SCIPgetNConss(scip);
2374  conss = SCIPgetConss(scip);
2375 
2376  /* initialize hash maps for keeping track of variables and constraints and their corresponding indices after being permuted by the ROC2-algorithm */
2377  SCIP_CALL( indexmapInit(detectordata->indexmap, vars, nvars, conss, nconss, detectordata->hashmapindices) );
2378 
2379 #ifdef WRITEALLOUTPUT
2380  {
2381  char filename[256];
2382  sprintf(filename, "%s_initial_problem", getProbNameWithoutPath(scip));
2383  //plotInitialProblem(scip, detectordata, filename);
2384  }
2385 #endif
2386 
2387  /* initialize index arrays ibegin, iend, jbegin, jend */
2388  SCIP_CALL( createRowindexList(scip, detectordata, detectordata->indexmap->indexcons, detectordata->indexmap->varindex, rowindices) );
2389 
2390  SCIP_CALL( createColumnindexList(scip, rowindices, columnindices) );
2391 
2392  SCIP_CALL( formIndexArray(detectordata->ibegin, detectordata->iend, rowindices) );
2393  SCIP_CALL( formIndexArray(detectordata->jbegin, detectordata->jend, columnindices) );
2394 
2395  /* ==================== */
2396  /* ===ROC2 algorithm=== */
2397  /* ==================== */
2398  SCIPdebugMessage("starting ROC2 algorithm\n");
2399 
2400 
2401 #ifdef WRITEALLOUTPUT
2402  SCIP_CALL( rankOrderClustering(scip, detectordata, detectordata->maxiterationsROC, &ROC_iterations) );
2403 
2404  /* check conditions for arrays ibegin and jbegin: ibegin[i]<=ibegin[i+k] for all positive k */
2405  if( ROC_iterations < detectordata->maxiterationsROC || detectordata->maxiterationsROC == -1 )
2406  {
2407  checkConsistencyOfIndexarrays(detectordata, nvars, nconss);
2408  }
2409 #else
2410  SCIP_CALL( rankOrderClustering(scip, detectordata, detectordata->maxiterationsROC, NULL) );
2411 #endif
2412  /* arrays jmin, jmax and minV */
2413  SCIPdebugMessage("calculating index arrays\n");
2414  detectordata->jmin[0] = detectordata->ibegin[0];
2415  detectordata->jmax[0] = detectordata->iend[0];
2416  detectordata->width[0] = detectordata->iend[0] - detectordata->ibegin[0];
2417  for( i = 1; i < nconss; ++i )
2418  {
2419  detectordata->width[i] = detectordata->iend[i] - detectordata->ibegin[i];
2420  detectordata->jmin[i] = detectordata->ibegin[i];
2421  detectordata->jmax[i] = MAX(detectordata->iend[i], detectordata->jmax[i-1]);
2422  detectordata->minV[i-1]=1 + (detectordata->jmax[i-1] - detectordata->jmin[i]);
2423  }
2424  /* ==================== */
2425  /* =====BLOCKING======= */
2426  /* ==================== */
2427 
2428  SCIP_CALL( blocking(scip, detectordata, decdecomps, ndecdecomps, nvars, nconss, result) );
2429  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, " found %d decompositions.\n", *ndecdecomps);
2430  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, " \tBlocks:", *ndecdecomps);
2431 #ifdef WRITEALLOUTPUT
2432  {
2433  char filename[256];
2434  sprintf(filename, "%s_ROC", getProbNameWithoutPath(scip));
2435  plotInitialProblem(scip, detectordata, filename);
2436  }
2437 #endif
2438  for( i = 0; i < *ndecdecomps; ++i )
2439  {
2440  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, " %i", DECdecompGetNBlocks( (*decdecomps)[i] ));
2441  }
2442  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, "\n");
2443 
2444  SCIP_CALL( SCIPreallocMemoryArray(scip, decdecomps, *ndecdecomps) );
2445 
2446  SCIP_CALL( freeData(scip, detectordata) );
2447 
2448  *result = SCIP_SUCCESS;
2449  return SCIP_OKAY;
2450 }
2451 
2452 static DEC_DECL_PROPAGATESEEED(detectorPropagateSeeedStairheur)
2453 {
2454  int i;
2455  int nconss; /* number of constraints in the problem */
2456  int nSeeeds;
2457  vector<vector<int> > rowindices;
2458  vector<vector<int> > columnindices;
2459  DEC_DETECTORDATA* detectordata = DECdetectorGetData(detector);
2460  gcg::Seeed* seeed;
2461 
2462  SCIP_CLOCK* temporaryClock;
2463  SCIP_CALL_ABORT(SCIPcreateClock(scip, &temporaryClock) );
2464  SCIP_CALL_ABORT( SCIPstartClock(scip, temporaryClock) );
2465 
2466  seeed = new gcg::Seeed(seeedPropagationData->seeedToPropagate);
2467  seeed->refineToMaster();
2468 
2469 #ifdef WRITEALLOUTPUT
2470  int ROC_iterations;
2471  SCIPwarningMessage(scip, "WRITEALLOUTPUT in detector stairheur is not implemented for seeeds.\n");
2472 #endif
2473 
2474  if( seeed->getNOpenconss() == 0 || seeed->getNOpenvars() == 0 )
2475  {
2476  delete seeed;
2477  seeedPropagationData->nNewSeeeds = 0;
2478  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock ) );
2479  SCIP_CALL_ABORT(SCIPfreeClock(scip, &temporaryClock) );
2480  *result = SCIP_SUCCESS;
2481  return SCIP_OKAY;
2482  }
2483 
2484  assert(scip != NULL);
2485  assert(detectordata != NULL);
2486 
2487  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, "Detecting stairheur structure:");
2488 
2489  SCIP_CALL( initData(scip, seeed, detectordata) );
2490 
2491  checkParameterConsistency(detectordata, result);
2492  nSeeeds = calculateNdecompositions(detectordata);
2493  SCIPdebugMessage("%i decompositions will be created\n", nSeeeds);
2494  seeedPropagationData->nNewSeeeds = 0;
2495 
2496  /* allocate space for output data */
2497  SCIP_CALL( SCIPallocMemoryArray(scip, &(seeedPropagationData->newSeeeds), nSeeeds) );
2498 
2499  nconss = seeed->getNOpenconss();
2500 
2501  /* initialize hash maps for keeping track of variables and constraints and their corresponding indices after being permuted by the ROC2-algorithm */
2502  SCIP_CALL( indexmapInit(detectordata->indexmap, seeed, detectordata->hashmapindices) );
2503 
2504 #ifdef WRITEALLOUTPUT
2505  {
2506  //char filename[256];
2507  //sprintf(filename, "%s_initial_problem", getProbNameWithoutPath(scip));
2508  }
2509 #endif
2510 
2511  /* initialize index arrays ibegin, iend, jbegin, jend */
2512  SCIP_CALL( createRowindexList(scip, seeed, seeedPropagationData->seeedpool, detectordata, detectordata->indexmap->indexcons, detectordata->indexmap->varindex, rowindices) );
2513 
2514  SCIP_CALL( createColumnindexList(scip, seeed, rowindices, columnindices) );
2515 
2516  SCIP_CALL( formIndexArray(detectordata->ibegin, detectordata->iend, rowindices) );
2517  SCIP_CALL( formIndexArray(detectordata->jbegin, detectordata->jend, columnindices) );
2518 
2519  /* ==================== */
2520  /* ===ROC2 algorithm=== */
2521  /* ==================== */
2522  SCIPdebugMessage("starting ROC2 algorithm\n");
2523 
2524 
2525 #ifdef WRITEALLOUTPUT
2526  SCIP_CALL( rankOrderClustering(scip, seeed, seeedPropagationData->seeedpool, detectordata, detectordata->maxiterationsROC, &ROC_iterations) );
2527 
2528  /* check conditions for arrays ibegin and jbegin: ibegin[i]<=ibegin[i+k] for all positive k */
2529  //if( ROC_iterations < detectordata->maxiterationsROC || detectordata->maxiterationsROC == -1 )
2530  //{
2531  // checkConsistencyOfIndexarrays(detectordata, seeed->getNOpenvars(), nconss);
2532  //}
2533 #else
2534  SCIP_CALL( rankOrderClustering(scip, seeed, seeedPropagationData->seeedpool, detectordata, detectordata->maxiterationsROC, NULL) );
2535 #endif
2536  /* arrays jmin, jmax and minV */
2537  SCIPdebugMessage("calculating index arrays\n");
2538  detectordata->jmin[0] = detectordata->ibegin[0];
2539  detectordata->jmax[0] = detectordata->iend[0];
2540  detectordata->width[0] = detectordata->iend[0] - detectordata->ibegin[0];
2541  for( i = 1; i < nconss; ++i )
2542  {
2543  detectordata->width[i] = detectordata->iend[i] - detectordata->ibegin[i];
2544  detectordata->jmin[i] = detectordata->ibegin[i];
2545  detectordata->jmax[i] = MAX(detectordata->iend[i], detectordata->jmax[i-1]);
2546  detectordata->minV[i-1]=1 + (detectordata->jmax[i-1] - detectordata->jmin[i]);
2547  }
2548  /* ==================== */
2549  /* =====BLOCKING======= */
2550  /* ==================== */
2551 
2552  SCIP_CALL_ABORT( SCIPstopClock(scip, temporaryClock ) );
2553  SCIP_CALL( blocking(scip, detectordata, seeed, seeedPropagationData->seeedpool, &(seeedPropagationData->newSeeeds), &(seeedPropagationData->nNewSeeeds), SCIPclockGetTime(temporaryClock), result) );
2554  SCIP_CALL_ABORT(SCIPfreeClock(scip, &temporaryClock) );
2555  SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, " found %d seeeds.\n", seeedPropagationData->nNewSeeeds);
2556  #ifdef WRITEALLOUTPUT
2557  {
2558  //char filename[256];
2559  //sprintf(filename, "%s_ROC", getProbNameWithoutPath(scip));
2560  //plotInitialProblem(scip, detectordata, filename);
2561  }
2562 #endif
2563  for( i = 0; i < seeedPropagationData->nNewSeeeds; ++i )
2564  {
2565  // SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, " %i", seeedPropagationData->newSeeeds[i]->getNBlocks());
2566  seeedPropagationData->newSeeeds[i]->setDetectorPropagated(detector);
2567  }
2568  // SCIPverbMessage(scip, SCIP_VERBLEVEL_NORMAL, NULL, "\n");
2569 
2570  SCIP_CALL( SCIPreallocMemoryArray(scip, &(seeedPropagationData->newSeeeds), seeedPropagationData->nNewSeeeds) );
2571 
2572 
2573  SCIP_CALL( freeData(scip, detectordata) );
2574 
2575  delete seeed;
2576  *result = SCIP_SUCCESS;
2577  return SCIP_OKAY;
2578 }
2579 #define detectorExitStairheur NULL
2580 #define detectorFinishSeeedStairheur NULL
2581 #define detectorPostprocessSeeedStairheur NULL
2582 
2583 #define setParamAggressiveStairheur NULL
2584 #define setParamDefaultStairheur NULL
2585 #define setParamFastStairheur NULL
2586 
2587 
2590 extern "C"
2592  SCIP* scip
2593  )
2594 {
2595  DEC_DETECTORDATA *detectordata = NULL;
2596  assert(scip != NULL);
2597 
2598  SCIP_CALL( SCIPallocMemory(scip, &detectordata) );
2599  assert(detectordata != NULL);
2600 
2601  detectordata->constoblock = NULL;
2602 
2604  detectordata, detectorDetectStairheur, detectorFreeStairheur, detectorInitStairheur, detectorExitStairheur, detectorPropagateSeeedStairheur, NULL, NULL, detectorFinishSeeedStairheur, detectorPostprocessSeeedStairheur, setParamAggressiveStairheur, setParamDefaultStairheur, setParamFastStairheur) );
2605 
2606 
2607 
2608  /* add stairheur detector parameters */
2609  SCIP_CALL( SCIPaddIntParam(scip, "detection/detectors/stairheur/nconssperblock",
2610  "The number of constraints per block (static blocking only)",
2611  &detectordata->nconssperblock, FALSE, DEFAULT_NCONSSPERBLOCK, 2, 1000000, NULL, NULL) );
2612  SCIP_CALL( SCIPaddIntParam(scip, "detection/detectors/stairheur/maxblocks",
2613  "The maximal number of blocks",
2614  &detectordata->maxblocks, FALSE, DEFAULT_MAXBLOCKS, 2, 1000000, NULL, NULL) );
2615  SCIP_CALL( SCIPaddIntParam(scip, "detection/detectors/stairheur/minblocks", "The minimal number of blocks",
2616  &detectordata->minblocks, FALSE, DEFAULT_MINBLOCKS, 2, 1000000, NULL, NULL) );
2617  SCIP_CALL( SCIPaddIntParam(scip, "detection/detectors/stairheur/desiredblocks",
2618  "The desired number of blocks. 0 means automatic determination of the number of blocks.",
2619  &detectordata->desiredblocks, FALSE, DEFAULT_DESIREDBLOCKS, 0, 1000000, NULL, NULL) );
2620  SCIP_CALL( SCIPaddBoolParam(scip, "detection/detectors/stairheur/dynamicblocking",
2621  "Enable blocking type 'dynamic'",
2622  &detectordata->dynamicblocking, FALSE, DEFAULT_DYNAMICBLOCKING, NULL, NULL) );
2623  SCIP_CALL( SCIPaddBoolParam(scip, "detection/detectors/stairheur/staticblocking",
2624  "Enable blocking type 'static'",
2625  &detectordata->staticblocking, FALSE, DEFAULT_STATICBLOCKING, NULL, NULL) );
2626  SCIP_CALL( SCIPaddBoolParam(scip, "detection/detectors/stairheur/blockingassoonaspossible",
2627  "Enable blocking type 'as soon as possible", &detectordata->blockingassoonaspossible,
2628  FALSE, DEFAULT_BLOCKINGASSOONASPOSSIBLE, NULL, NULL) );
2629  SCIP_CALL( SCIPaddBoolParam(scip, "detection/detectors/stairheur/multipledecomps",
2630  "Enables multiple decompositions for all enabled blocking types. Ranging from minblocks to maxblocks",
2631  &detectordata->multipledecomps, FALSE, DEFAULT_MULTIPLEDECOMPS, NULL, NULL) );
2632  SCIP_CALL( SCIPaddIntParam(scip, "detection/detectors/stairheur/maxiterationsROC",
2633  "The maximum number of iterations of the ROC-algorithm. -1 for no limit",
2634  &detectordata->maxiterationsROC, FALSE, DEFAULT_MAXITERATIONSROC, -1, 1000000, NULL, NULL) );
2635 
2636  return SCIP_OKAY;
2637 }
#define DEC_ENABLEDPOSTPROCESSING
static SCIP_RETCODE resetDetectordata(SCIP *scip, DEC_DETECTORDATA *detectordata)
vector< int > * blockedAfterrow
static DEC_DECL_FREEDETECTOR(detectorFreeStairheur)
SCIP_Bool staticblocking
SCIP_Bool dynamicblocking
static SCIP_RETCODE indexmapCreate(SCIP *scip, INDEXMAP **indexmap, int nconss, int nvars)
int getNVarsForCons(int consIndex)
returns the number of variables for a given constraint
SCIP_RETCODE SCIPincludeDetectorStairheur(SCIP *scip)
const int * getOpenconss()
returns array containing constraints not assigned yet
static void checkParameterConsistency(DEC_DETECTORDATA *detectordata, SCIP_RESULT *result)
static vector< int >::iterator nextRowToBlockAt(DEC_DETECTORDATA *detectordata, vector< int >::iterator it_constrictions, vector< int > *it_vector, int min_block_size, int prev_block_first_row, int prev_block_last_row)
static vector< int >::iterator findBlockingCandidate(vector< int >::iterator it_constrictions, vector< int > *it_vector, int min_block_size, int prev_block_last_row)
GCG interface methods.
DEC_DETECTORDATA * DECdetectorGetData(DEC_DETECTOR *detector)
returns the data of the provided detector
static SCIP_RETCODE formIndexArray(int *begin, int *end, vector< vector< int > > &indices)
static int getMinColIndex(DEC_DETECTORDATA *detectordata, int row)
#define detectorPostprocessSeeedStairheur
static int getMaxColIndex(DEC_DETECTORDATA *detectordata, int from_row, int to_row)
static void SCIPvectorRearrange(vector< vector< int > > &l, vector< int > order)
static SCIP_RETCODE initData(SCIP *scip, DEC_DETECTORDATA *detectordata)
int GCGconsGetNVars(SCIP *scip, SCIP_CONS *cons)
Definition: scip_misc.c:433
static SCIP_RETCODE rowsWithConstriction(SCIP *scip, DEC_DETECTORDATA *detectordata)
#define DEC_LEGACYMODE
int getNOpenvars()
returns size of vector containing variables not assigned yet
#define DEFAULT_MULTIPLEDECOMPS
void printvector(vector< int > l)
static DEC_DECL_INITDETECTOR(detectorInitStairheur)
SCIP_RETCODE refineToMaster()
refine seeed with focus on master: do obvious (
#define DEFAULT_MINBLOCKS
#define DEC_ENABLEDFINISHING
static SCIP_RETCODE rankOrderClustering(SCIP *scip, DEC_DETECTORDATA *detectordata, int max_iterations, int *iterations)
static SCIP_RETCODE freeData(SCIP *scip, DEC_DETECTORDATA *detectordata)
detector for staircase structures via ROC algorithms
static SCIP_RETCODE rankOrderClusteringIteration(SCIP *scip, DEC_DETECTORDATA *detectordata, INDEXMAP *inputmap, INDEXMAP *outputmap)
static SCIP_RETCODE blockingAsSoonAsPossible(SCIP *scip, DEC_DETECTORDATA *detectordata, int desired_blocks, int nvars)
#define setParamAggressiveStairheur
bool isVarOpenvar(int var)
returns true if the var is an open var
static SCIP_RETCODE blockingStatic(SCIP *scip, DEC_DETECTORDATA *detectordata)
#define DEFAULT_DESIREDBLOCKS
static DEC_DECL_PROPAGATESEEED(detectorPropagateSeeedStairheur)
static SCIP_Bool isValidBlocking(DEC_DETECTORDATA *detectordata, int prev_block_first_row, int prev_block_last_row, int block_at_row)
int DECdecompGetNBlocks(DEC_DECOMP *decomp)
Definition: decomp.c:740
static SCIP_RETCODE blockingDynamic(SCIP *scip, DEC_DETECTORDATA *detectordata, int tau, int nvars)
static SCIP_RETCODE blocking(SCIP *scip, DEC_DETECTORDATA *detectordata, DEC_DECOMP ***decdecomps, int *ndecdecomps, int nvars, int ncons, SCIP_RESULT *result)
#define DEC_FREQCALLROUNDORIGINAL
various SCIP helper methods
class to manage partial decompositions (aka seeed), each seeed corresponds to one seeedpool which con...
Definition: class_seeed.h:71
#define DEFAULT_BLOCKINGASSOONASPOSSIBLE
SCIP_HASHMAP * varindex
static DEC_DECL_DETECTSTRUCTURE(detectorDetectStairheur)
#define DEFAULT_NCONSSPERBLOCK
SCIP_RETCODE DECdecompCreate(SCIP *scip, DEC_DECOMP **decdecomp)
Definition: decomp.c:467
SCIP_HASHMAP * indexvar
static SCIP_Bool arraysAreEqual(int *new_array, int *old_array, int num_elements)
SCIP_RETCODE GCGconsGetVars(SCIP *scip, SCIP_CONS *cons, SCIP_VAR **vars, int nvars)
Definition: scip_misc.c:489
#define DEC_SKIP
static SCIP_RETCODE createColumnindexList(SCIP *scip, vector< vector< int > > &rowindices, vector< vector< int > > &columnindices)
SCIP_HASHMAP * consindex
#define DEC_MINCALLROUNDORIGINAL
#define DEC_DECCHAR
static void indexmapFree(SCIP *scip, INDEXMAP **indexmap)
#define DEC_DESC
SCIP_Bool multipledecomps
SCIP_HASHMAP * constoblock
#define detectorExitStairheur
#define DEC_MINCALLROUND
static vector< int > SCIPvectorCreateInt(int from, int to)
#define DEFAULT_DYNAMICBLOCKING
vector< int > * rowsWithConstrictions
#define DEC_ENABLEDORIGINAL
#define DEC_MAXCALLROUND
#define setParamDefaultStairheur
#define DEC_PRIORITY
static SCIP_RETCODE createRowindexList(SCIP *scip, DEC_DETECTORDATA *detectordata, SCIP_HASHMAP *indexcons, SCIP_HASHMAP *varindex, vector< vector< int > > &rowindices)
static SCIP_RETCODE assignConsToBlock(SCIP *scip, DEC_DETECTORDATA *detectordata, int block, int first_cons, int last_cons)
SCIP_RETCODE DECincludeDetector(SCIP *scip, const char *name, const char decchar, const char *description, int freqCallRound, int maxCallRound, int minCallRound, int freqCallRoundOriginal, int maxCallRoundOriginal, int minCallRoundOriginal, int priority, SCIP_Bool enabled, SCIP_Bool enabledOriginal, SCIP_Bool enabledFinishing, SCIP_Bool enabledPostprocessing, SCIP_Bool skip, SCIP_Bool usefulRecall, SCIP_Bool legacymode, DEC_DETECTORDATA *detectordata, DEC_DECL_DETECTSTRUCTURE((*detectStructure)), DEC_DECL_FREEDETECTOR((*freeDetector)), DEC_DECL_INITDETECTOR((*initDetector)), DEC_DECL_EXITDETECTOR((*exitDetector)), DEC_DECL_PROPAGATESEEED((*propagateSeeedDetector)), DEC_DECL_PROPAGATEFROMTOOLBOX((*propagateFromToolboxDetector)), DEC_DECL_FINISHFROMTOOLBOX((*finishFromToolboxDetector)), DEC_DECL_FINISHSEEED((*finishSeeedDetector)), DEC_DECL_POSTPROCESSSEEED((*postprocessSeeedDetector)), DEC_DECL_SETPARAMAGGRESSIVE((*setParamAggressiveDetector)), DEC_DECL_SETPARAMDEFAULT((*setParamDefaultDetector)),)
includes one detector
static vector< int > rowOrdering(SCIP *scip, vector< vector< int > > &columnindices, int nrows)
SCIP_Bool blockingassoonaspossible
const int * getOpenvars()
returns array containing variables not assigned yet
#define DEFAULT_STATICBLOCKING
#define DEC_FREQCALLROUND
#define DEFAULT_MAXBLOCKS
#define DEC_USEFULRECALL
static SCIP_RETCODE indexmapInit(INDEXMAP *indexmap, SCIP_VAR **vars, int nvars, SCIP_CONS **conss, int nconss, int *hashmapindices)
SCIP_RETCODE DECfilloutDecompFromConstoblock(SCIP *scip, DEC_DECOMP *decomp, SCIP_HASHMAP *constoblock, int nblocks, SCIP_Bool staircase)
Definition: decomp.c:1454
#define DEC_MAXCALLROUNDORIGINAL
#define DEC_DETECTORNAME
#define setParamFastStairheur
int getNOpenconss()
returns size of vector containing constraints not assigned yet
class with functions for seeed pool where a seeed is a (potentially incomplete) description of a deco...
const char * DECdetectorGetName(DEC_DETECTOR *detector)
returns the name of the provided detector
#define DEFAULT_MAXITERATIONSROC
#define detectorFinishSeeedStairheur
constraint handler for structure detection
struct DecDecomp DEC_DECOMP
Definition: type_decomp.h:43
#define DEC_ENABLED
static int calculateNdecompositions(DEC_DETECTORDATA *detectordata)
public methods for working with decomposition structures
void printnested(vector< vector< int > > l)
const int * getVarsForCons(int consIndex)
returns the variable indices of the coefficient matrix for a constraint
SCIP_HASHMAP * indexcons