# SCIP

Solving Constraint Integer Programs

cons_cumulative.c
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2019 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scip.zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15
16 /**@file cons_cumulative.c
17  * @brief constraint handler for cumulative constraints
18  * @author Timo Berthold
19  * @author Stefan Heinz
20  * @author Jens Schulz
21  *
22  * Given:
23  * - a set of jobs, represented by their integer start time variables \f$S_j\f$, their array of processing times \f$p_j\f$ and of
24  * their demands \f$d_j\f$.
25  * - an integer resource capacity \f$C\f$
26  *
27  * The cumulative constraint ensures that for each point in time \f$t\f$ \f$\sum_{j: S_j \leq t < S_j + p_j} d_j \leq C\f$ holds.
28  *
29  * Separation:
30  * - can be done using binary start time model, see Pritskers, Watters and Wolfe
31  * - or by just separating relatively weak cuts on the integer start time variables
32  *
33  * Propagation:
34  * - time tabling, Klein & Scholl (1999)
35  * - Edge-finding from Petr Vilim, adjusted and simplified for dynamic repropagation
36  * (2009)
37  * - energetic reasoning, see Baptiste, Le Pape, Nuijten (2001)
38  *
39  */
40
41 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
42
43 #include <assert.h>
44 #include <string.h>
45
46 #include "tclique/tclique.h"
47 #include "scip/cons_cumulative.h"
49 #include "scip/cons_knapsack.h"
50 #include "scip/scipdefplugins.h"
51
52 /**@name Constraint handler properties
53  *
54  * @{
55  */
56
57 /* constraint handler properties */
58 #define CONSHDLR_NAME "cumulative"
59 #define CONSHDLR_DESC "cumulative constraint handler"
60 #define CONSHDLR_SEPAPRIORITY 2100000 /**< priority of the constraint handler for separation */
61 #define CONSHDLR_ENFOPRIORITY -2040000 /**< priority of the constraint handler for constraint enforcing */
62 #define CONSHDLR_CHECKPRIORITY -3030000 /**< priority of the constraint handler for checking feasibility */
63 #define CONSHDLR_SEPAFREQ 1 /**< frequency for separating cuts; zero means to separate only in the root node */
64 #define CONSHDLR_PROPFREQ 1 /**< frequency for propagating domains; zero means only preprocessing propagation */
65 #define CONSHDLR_EAGERFREQ 100 /**< frequency for using all instead of only the useful constraints in separation,
66  * propagation and enforcement, -1 for no eager evaluations, 0 for first only */
67 #define CONSHDLR_MAXPREROUNDS -1 /**< maximal number of presolving rounds the constraint handler participates in (-1: no limit) */
68 #define CONSHDLR_DELAYSEPA FALSE /**< should separation method be delayed, if other separators found cuts? */
69 #define CONSHDLR_DELAYPROP FALSE /**< should propagation method be delayed, if other propagators found reductions? */
70 #define CONSHDLR_NEEDSCONS TRUE /**< should the constraint handler be skipped, if no constraints are available? */
72 #define CONSHDLR_PRESOLTIMING SCIP_PRESOLTIMING_ALWAYS
73 #define CONSHDLR_PROP_TIMING SCIP_PROPTIMING_BEFORELP
75 /**@} */
76
77 /**@name Default parameter values
78  *
79  * @{
80  */
81
82 /* default parameter values */
83
84 /* separation */
85 #define DEFAULT_USEBINVARS FALSE /**< should the binary representation be used? */
86 #define DEFAULT_LOCALCUTS FALSE /**< should cuts be added only locally? */
87 #define DEFAULT_USECOVERCUTS TRUE /**< should covering cuts be added? */
88 #define DEFAULT_CUTSASCONSS TRUE /**< should the cuts be created as knapsack constraints? */
89 #define DEFAULT_SEPAOLD TRUE /**< shall old sepa algo be applied? */
91 /* propagation */
92 #define DEFAULT_TTINFER TRUE /**< should time-table (core-times) propagator be used to infer bounds? */
93 #define DEFAULT_EFCHECK FALSE /**< should edge-finding be used to detect an overload? */
94 #define DEFAULT_EFINFER FALSE /**< should edge-finding be used to infer bounds? */
95 #define DEFAULT_USEADJUSTEDJOBS FALSE /**< should during edge-finding jobs be adusted which run on the border of the effective time horizon? */
96 #define DEFAULT_TTEFCHECK TRUE /**< should time-table edge-finding be used to detect an overload? */
97 #define DEFAULT_TTEFINFER TRUE /**< should time-table edge-finding be used to infer bounds? */
99 /* presolving */
100 #define DEFAULT_DUALPRESOLVE TRUE /**< should dual presolving be applied? */
101 #define DEFAULT_COEFTIGHTENING FALSE /**< should coeffisient tightening be applied? */
102 #define DEFAULT_NORMALIZE TRUE /**< should demands and capacity be normalized? */
103 #define DEFAULT_PRESOLPAIRWISE TRUE /**< should pairwise constraint comparison be performed in presolving? */
104 #define DEFAULT_DISJUNCTIVE TRUE /**< extract disjunctive constraints? */
105 #define DEFAULT_DETECTDISJUNCTIVE TRUE /**< search for conflict set via maximal cliques to detect disjunctive constraints */
106 #define DEFAULT_DETECTVARBOUNDS TRUE /**< search for conflict set via maximal cliques to detect variable bound constraints */
107 #define DEFAULT_MAXNODES 10000LL /**< number of branch-and-bound nodes to solve an independent cumulative constraint (-1: no limit) */
109 /* enforcement */
110 #define DEFAULT_FILLBRANCHCANDS FALSE /**< should branching candidates be added to storage? */
112 /* conflict analysis */
113 #define DEFAULT_USEBDWIDENING TRUE /**< should bound widening be used during conflict analysis? */
115 /**@} */
116
117 /**@name Event handler properties
118  *
119  * @{
120  */
121
122 #define EVENTHDLR_NAME "cumulative"
123 #define EVENTHDLR_DESC "bound change event handler for cumulative constraints"
125 /**@} */
126
127 /*
128  * Data structures
129  */
130
131 /** constraint data for cumulative constraints */
132 struct SCIP_ConsData
133 {
134  SCIP_VAR** vars; /**< array of variable representing the start time of each job */
135  SCIP_Bool* downlocks; /**< array to store if the variable has a down lock */
136  SCIP_Bool* uplocks; /**< array to store if the variable has an uplock */
137  SCIP_CONS** linkingconss; /**< array of linking constraints for the integer variables */
138  SCIP_ROW** demandrows; /**< array of rows of linear relaxation of this problem */
139  SCIP_ROW** scoverrows; /**< array of rows of small cover cuts of this problem */
140  SCIP_ROW** bcoverrows; /**< array of rows of big cover cuts of this problem */
141  int* demands; /**< array containing corresponding demands */
142  int* durations; /**< array containing corresponding durations */
143  SCIP_Real resstrength1; /**< stores the resource strength 1*/
144  SCIP_Real resstrength2; /**< stores the resource strength 2 */
145  SCIP_Real cumfactor1; /**< stroes the cumulativeness of the constraint */
146  SCIP_Real disjfactor1; /**< stores the disjunctiveness of the constraint */
147  SCIP_Real disjfactor2; /**< stores the disjunctiveness of the constraint */
148  SCIP_Real estimatedstrength;
149  int nvars; /**< number of variables */
150  int varssize; /**< size of the arrays */
151  int ndemandrows; /**< number of rows of cumulative constrint for linear relaxation */
152  int demandrowssize; /**< size of array rows of demand rows */
153  int nscoverrows; /**< number of rows of small cover cuts */
154  int scoverrowssize; /**< size of array of small cover cuts */
155  int nbcoverrows; /**< number of rows of big cover cuts */
156  int bcoverrowssize; /**< size of array of big cover cuts */
157  int capacity; /**< available cumulative capacity */
158
159  int hmin; /**< left bound of time axis to be considered (including hmin) */
160  int hmax; /**< right bound of time axis to be considered (not including hmax) */
161
162  unsigned int signature; /**< constraint signature which is need for pairwise comparison */
163
164  unsigned int validsignature:1; /**< is the signature valid */
165  unsigned int normalized:1; /**< is the constraint normalized */
166  unsigned int covercuts:1; /**< cover cuts are created? */
167  unsigned int propagated:1; /**< is constraint propagted */
168  unsigned int varbounds:1; /**< bool to store if variable bound strengthening was already preformed */
169  unsigned int triedsolving:1; /**< bool to store if we tried already to solve that constraint as independent subproblem */
170
171 #ifdef SCIP_STATISTIC
172  int maxpeak;
173 #endif
174 };
175
176 /** constraint handler data */
177 struct SCIP_ConshdlrData
178 {
179  SCIP_EVENTHDLR* eventhdlr; /**< event handler for bound change events */
180
181  SCIP_Bool usebinvars; /**< should the binary variables be used? */
182  SCIP_Bool cutsasconss; /**< should the cumulative constraint create cuts as knapsack constraints? */
183  SCIP_Bool ttinfer; /**< should time-table (core-times) propagator be used to infer bounds? */
184  SCIP_Bool efcheck; /**< should edge-finding be used to detect an overload? */
185  SCIP_Bool efinfer; /**< should edge-finding be used to infer bounds? */
186  SCIP_Bool useadjustedjobs; /**< should during edge-finding jobs be adusted which run on the border of the effective time horizon? */
187  SCIP_Bool ttefcheck; /**< should time-table edge-finding be used to detect an overload? */
188  SCIP_Bool ttefinfer; /**< should time-table edge-finding be used to infer bounds? */
189  SCIP_Bool localcuts; /**< should cuts be added only locally? */
190  SCIP_Bool usecovercuts; /**< should covering cuts be added? */
191  SCIP_Bool sepaold; /**< shall old sepa algo be applied? */
192
193  SCIP_Bool fillbranchcands; /**< should branching candidates be added to storage? */
194
195  SCIP_Bool dualpresolve; /**< should dual presolving be applied? */
196  SCIP_Bool coeftightening; /**< should coeffisient tightening be applied? */
197  SCIP_Bool normalize; /**< should demands and capacity be normalized? */
198  SCIP_Bool disjunctive; /**< extract disjunctive constraints? */
199  SCIP_Bool detectdisjunctive; /**< search for conflict set via maximal cliques to detect disjunctive constraints */
200  SCIP_Bool detectvarbounds; /**< search for conflict set via maximal cliques to detect variable bound constraints */
201  SCIP_Bool usebdwidening; /**< should bound widening be used during conflict analysis? */
202  SCIP_Bool presolpairwise; /**< should pairwise constraint comparison be performed in presolving? */
203  SCIP_Bool detectedredundant; /**< was detection of redundant constraints already performed? */
204
205  SCIP_Longint maxnodes; /**< number of branch-and-bound nodes to solve an independent cumulative constraint (-1: no limit) */
206
207  SCIP_DECL_SOLVECUMULATIVE((*solveCumulative)); /**< method to use a single cumulative condition */
208
209  /* statistic values which are collected if SCIP_STATISTIC is defined */
210 #ifdef SCIP_STATISTIC
211  SCIP_Longint nlbtimetable; /**< number of times the lower bound was tightened by the time-table propagator */
212  SCIP_Longint nubtimetable; /**< number of times the upper bound was tightened by the time-table propagator */
213  SCIP_Longint ncutofftimetable; /**< number of times the a cutoff was detected due to time-table propagator */
214  SCIP_Longint nlbedgefinder; /**< number of times the lower bound was tightened by the edge-finder propagator */
215  SCIP_Longint nubedgefinder; /**< number of times the upper bound was tightened by the edge-finder propagator */
216  SCIP_Longint ncutoffedgefinder; /**< number of times the a cutoff was detected due to edge-finder propagator */
217  SCIP_Longint ncutoffoverload; /**< number of times the a cutoff was detected due to overload checking via edge-finding */
218  SCIP_Longint nlbTTEF; /**< number of times the lower bound was tightened by time-table edge-finding */
219  SCIP_Longint nubTTEF; /**< number of times the upper bound was tightened by time-table edge-finding */
220  SCIP_Longint ncutoffoverloadTTEF;/**< number of times the a cutoff was detected due to overload checking via time-table edge-finding */
221
222  int nirrelevantjobs; /**< number of time a irrelevant/redundant jobs was removed form a constraint */
223  int nalwaysruns; /**< number of time a job removed form a constraint which run completely during the effective horizon */
224  int nremovedlocks; /**< number of times a up or down lock was removed */
225  int ndualfixs; /**< number of times a dual fix was performed by a single constraint */
226  int ndecomps; /**< number of times a constraint was decomposed */
227  int ndualbranchs; /**< number of times a dual branch was discoverd and applicable via probing */
228  int nallconsdualfixs; /**< number of times a dual fix was performed due to knowledge of all cumulative constraints */
229  int naddedvarbounds; /**< number of added variable bounds constraints */
230  int naddeddisjunctives; /**< number of added disjunctive constraints */
231
232  SCIP_Bool iscopy; /**< Boolean to store if constraint handler is part of a copy */
233 #endif
234 };
235
236 /**@name Inference Information Methods
237  *
238  * An inference information can be passed with each domain reduction to SCIP. This information is passed back to the
239  * constraint handler if the corresponding bound change has to be explained. It can be used to store information which
240  * help to construct a reason/explanation for a bound change. The inference information is limited to size of integer.
241  *
242  * In case of the cumulative constraint handler we store the used propagation algorithms for that particular bound
243  * change and the earliest start and latest completion time of all jobs in the conflict set.
244  *
245  * @{
246  */
247
248 /** Propagation rules */
249 enum Proprule
250 {
251  PROPRULE_1_CORETIMES = 1, /**< core-time propagator */
252  PROPRULE_2_EDGEFINDING = 2, /**< edge-finder */
253  PROPRULE_3_TTEF = 3 /**< time-table edeg-finding */
254 };
255 typedef enum Proprule PROPRULE;
257 /** inference information */
258 struct InferInfo
259 {
260  union
261  {
262  /** struct to use the inference information */
263  struct
264  {
265  unsigned int proprule:2; /**< propagation rule that was applied */
266  unsigned int data1:15; /**< data field one */
267  unsigned int data2:15; /**< data field two */
268  } asbits;
269  int asint; /**< inference information as a single int value */
270  } val;
271 };
272 typedef struct InferInfo INFERINFO;
274 /** converts an integer into an inference information */
275 static
276 INFERINFO intToInferInfo(
277  int i /**< integer to convert */
278  )
279 {
280  INFERINFO inferinfo;
281
282  inferinfo.val.asint = i;
283
284  return inferinfo;
285 }
286
287 /** converts an inference information into an int */
288 static
289 int inferInfoToInt(
290  INFERINFO inferinfo /**< inference information to convert */
291  )
292 {
293  return inferinfo.val.asint;
294 }
295
296 /** returns the propagation rule stored in the inference information */
297 static
299  INFERINFO inferinfo /**< inference information to convert */
300  )
301 {
302  return (PROPRULE) inferinfo.val.asbits.proprule;
303 }
304
305 /** returns data field one of the inference information */
306 static
308  INFERINFO inferinfo /**< inference information to convert */
309  )
310 {
311  return (int) inferinfo.val.asbits.data1;
312 }
313
314 /** returns data field two of the inference information */
315 static
317  INFERINFO inferinfo /**< inference information to convert */
318  )
319 {
320  return (int) inferinfo.val.asbits.data2;
321 }
322
323
324 /** constructs an inference information out of a propagation rule, an earliest start and a latest completion time */
325 static
326 INFERINFO getInferInfo(
327  PROPRULE proprule, /**< propagation rule that deduced the value */
328  int data1, /**< data field one */
329  int data2 /**< data field two */
330  )
331 {
332  INFERINFO inferinfo;
333
334  /* check that the data menber are in the range of the available bits */
335  assert((int)proprule <= (1<<2)-1); /*lint !e685*/
336  assert(data1 >= 0 && data1 < (1<<15));
337  assert(data2 >= 0 && data2 < (1<<15));
338
339  inferinfo.val.asbits.proprule = proprule; /*lint !e641*/
340  inferinfo.val.asbits.data1 = (unsigned int) data1; /*lint !e732*/
341  inferinfo.val.asbits.data2 = (unsigned int) data2; /*lint !e732*/
342
343  return inferinfo;
344 }
345
346 /**@} */
347
348 /*
349  * Local methods
350  */
351
352 /**@name Miscellaneous Methods
353  *
354  * @{
355  */
356
357 #ifndef NDEBUG
358
359 /** compute the core of a job which lies in certain interval [begin, end) */
360 static
362  int begin, /**< begin of the interval */
363  int end, /**< end of the interval */
364  int ect, /**< earliest completion time */
365  int lst /**< latest start time */
366  )
367 {
368  int core;
369
370  core = MAX(0, MIN(end, ect) - MAX(lst, begin));
371
372  return core;
373 }
374 #else
375 #define computeCoreWithInterval(begin, end, ect, lst) (MAX(0, MIN((end), (ect)) - MAX((lst), (begin))))
376 #endif
377
378 /** returns the implied earliest start time */
379 static
381  SCIP* scip, /**< SCIP data structure */
382  SCIP_VAR* var, /**< variable for which the implied est should be returned */
383  SCIP_HASHMAP* addedvars, /**< hash map containig the variable which are already added */
384  int* est /**< pointer to store the implied earliest start time */
385  )
386 { /*lint --e{715}*/
387 #if 0
388  SCIP_VAR** vbdvars;
389  SCIP_VAR* vbdvar;
390  SCIP_Real* vbdcoefs;
391  SCIP_Real* vbdconsts;
392  void* image;
393  int nvbdvars;
394  int v;
395 #endif
396
397  (*est) = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
398
399 #if 0
400  /* the code contains a bug; we need to check if an implication forces that the jobs do not run in parallel */
401
402  nvbdvars = SCIPvarGetNVlbs(var);
403  vbdvars = SCIPvarGetVlbVars(var);
404  vbdcoefs = SCIPvarGetVlbCoefs(var);
405  vbdconsts = SCIPvarGetVlbConstants(var);
406
407  for( v = 0; v < nvbdvars; ++v )
408  {
409  vbdvar = vbdvars[v];
410  assert(vbdvar != NULL);
411
412  image = SCIPhashmapGetImage(addedvars, (void*)vbdvar);
413
414  if( image != NULL && SCIPisEQ(scip, vbdcoefs[v], 1.0 ) )
415  {
416  int duration;
417  int vbdconst;
418
419  duration = (int)(size_t)image;
420  vbdconst = SCIPconvertRealToInt(scip, vbdconsts[v]);
421
422  SCIPdebugMsg(scip, "check implication <%s>[%g,%g] >= <%s>[%g,%g] + <%g>\n",
424  SCIPvarGetName(vbdvar), SCIPvarGetLbLocal(vbdvar), SCIPvarGetUbLocal(vbdvar), vbdconsts[v]);
425
426  if( duration >= vbdconst )
427  {
428  int impliedest;
429
430  impliedest = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vbdvar)) + duration;
431
432  if( (*est) < impliedest )
433  {
434  (*est) = impliedest;
435
436  SCIP_CALL( SCIPhashmapRemove(addedvars, (void*)vbdvar) );
437  }
438  }
439  }
440  }
441 #endif
442
443  return SCIP_OKAY;
444 }
445
446 /** returns the implied latest completion time */
447 static
449  SCIP* scip, /**< SCIP data structure */
450  SCIP_VAR* var, /**< variable for which the implied est should be returned */
451  int duration, /**< duration of the given job */
452  SCIP_HASHMAP* addedvars, /**< hash map containig the variable which are already added */
453  int* lct /**< pointer to store the implied latest completion time */
454  )
455 { /*lint --e{715}*/
456 #if 0
457  SCIP_VAR** vbdvars;
458  SCIP_VAR* vbdvar;
459  SCIP_Real* vbdcoefs;
460  SCIP_Real* vbdconsts;
461  int nvbdvars;
462  int v;
463 #endif
464
465  (*lct) = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var)) + duration;
466
467 #if 0
468  /* the code contains a bug; we need to check if an implication forces that the jobs do not run in parallel */
469
470  nvbdvars = SCIPvarGetNVubs(var);
471  vbdvars = SCIPvarGetVubVars(var);
472  vbdcoefs = SCIPvarGetVubCoefs(var);
473  vbdconsts = SCIPvarGetVubConstants(var);
474
475  for( v = 0; v < nvbdvars; ++v )
476  {
477  vbdvar = vbdvars[v];
478  assert(vbdvar != NULL);
479
480  if( SCIPhashmapExists(addedvars, (void*)vbdvar) && SCIPisEQ(scip, vbdcoefs[v], 1.0 ) )
481  {
482  int vbdconst;
483
484  vbdconst = SCIPconvertRealToInt(scip, -vbdconsts[v]);
485
486  SCIPdebugMsg(scip, "check implication <%s>[%g,%g] <= <%s>[%g,%g] + <%g>\n",
488  SCIPvarGetName(vbdvar), SCIPvarGetLbLocal(vbdvar), SCIPvarGetUbLocal(vbdvar), vbdconsts[v]);
489
490  if( duration >= -vbdconst )
491  {
492  int impliedlct;
493
494  impliedlct = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vbdvar));
495
496  if( (*lct) > impliedlct )
497  {
498  (*lct) = impliedlct;
499
500  SCIP_CALL( SCIPhashmapRemove(addedvars, (void*)vbdvar) );
501  }
502  }
503  }
504  }
505 #endif
506
507  return SCIP_OKAY;
508 }
509
510 /** collects all necessary binary variables to represent the jobs which can be active at time point of interest */
511 static
513  SCIP* scip, /**< SCIP data structure */
514  SCIP_CONSDATA* consdata, /**< constraint data */
515  SCIP_VAR*** vars, /**< pointer to the array to store the binary variables */
516  int** coefs, /**< pointer to store the coefficients */
517  int* nvars, /**< number if collect binary variables */
518  int* startindices, /**< permutation with rspect to the start times */
519  int curtime, /**< current point in time */
520  int nstarted, /**< number of jobs that start before the curtime or at curtime */
521  int nfinished /**< number of jobs that finished before curtime or at curtime */
522  )
523 {
524  int nrowvars;
525  int startindex;
526  int size;
527
528  size = 10;
529  nrowvars = 0;
530  startindex = nstarted - 1;
531
532  SCIP_CALL( SCIPallocBufferArray(scip, vars, size) );
533  SCIP_CALL( SCIPallocBufferArray(scip, coefs, size) );
534
535  /* search for the (nstarted - nfinished) jobs which are active at curtime */
536  while( nstarted - nfinished > nrowvars )
537  {
538  SCIP_VAR* var;
539  int endtime;
540  int duration;
541  int demand;
542  int varidx;
543
544  /* collect job information */
545  varidx = startindices[startindex];
546  assert(varidx >= 0 && varidx < consdata->nvars);
547
548  var = consdata->vars[varidx];
549  duration = consdata->durations[varidx];
550  demand = consdata->demands[varidx];
551  assert(var != NULL);
552
553  endtime = SCIPconvertRealToInt(scip, SCIPvarGetUbGlobal(var)) + duration;
554
555  /* check the end time of this job is larger than the curtime; in this case the job is still running */
556  if( endtime > curtime )
557  {
558  SCIP_VAR** binvars;
559  int* vals;
560  int nbinvars;
561  int start;
562  int end;
563  int b;
564
565  /* check if the linking constraints exists */
567  assert(SCIPgetConsLinking(scip, var) != NULL);
569
570  /* collect linking constraint information */
571  SCIP_CALL( SCIPgetBinvarsLinking(scip, consdata->linkingconss[varidx], &binvars, &nbinvars) );
573
574  start = curtime - duration + 1;
575  end = MIN(curtime, endtime - duration);
576
577  for( b = 0; b < nbinvars; ++b )
578  {
579  if( vals[b] < start )
580  continue;
581
582  if( vals[b] > end )
583  break;
584
585  assert(binvars[b] != NULL);
586
587  /* ensure array proper array size */
588  if( size == *nvars )
589  {
590  size *= 2;
591  SCIP_CALL( SCIPreallocBufferArray(scip, vars, size) );
592  SCIP_CALL( SCIPreallocBufferArray(scip, coefs, size) );
593  }
594
595  (*vars)[*nvars] = binvars[b];
596  (*coefs)[*nvars] = demand;
597  (*nvars)++;
598  }
599  nrowvars++;
600  }
601
602  startindex--;
603  }
604
605  return SCIP_OKAY;
606 }
607
608 /** collect all integer variable which belong to jobs which can run at the point of interest */
609 static
611  SCIP* scip, /**< SCIP data structure */
612  SCIP_CONSDATA* consdata, /**< constraint data */
613  SCIP_VAR*** activevars, /**< jobs that are currently running */
614  int* startindices, /**< permutation with rspect to the start times */
615  int curtime, /**< current point in time */
616  int nstarted, /**< number of jobs that start before the curtime or at curtime */
617  int nfinished, /**< number of jobs that finished before curtime or at curtime */
618  SCIP_Bool lower, /**< shall cuts be created due to lower or upper bounds? */
619  int* lhs /**< lhs for the new row sum of lbs + minoffset */
620  )
621 {
622  SCIP_VAR* var;
623  int startindex;
624  int endtime;
625  int duration;
626  int starttime;
627
628  int varidx;
629  int sumofstarts;
630  int mindelta;
631  int counter;
632
633  assert(curtime >= consdata->hmin);
634  assert(curtime < consdata->hmax);
635
636  counter = 0;
637  sumofstarts = 0;
638
639  mindelta = INT_MAX;
640
641  startindex = nstarted - 1;
642
643  /* search for the (nstarted - nfinished) jobs which are active at curtime */
644  while( nstarted - nfinished > counter )
645  {
646  assert(startindex >= 0);
647
648  /* collect job information */
649  varidx = startindices[startindex];
650  assert(varidx >= 0 && varidx < consdata->nvars);
651
652  var = consdata->vars[varidx];
653  duration = consdata->durations[varidx];
654  assert(duration > 0);
655  assert(var != NULL);
656
657  if( lower )
658  starttime = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
659  else
660  starttime = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
661
662  endtime = MIN(starttime + duration, consdata->hmax);
663
664  /* check the end time of this job is larger than the curtime; in this case the job is still running */
665  if( endtime > curtime )
666  {
667  (*activevars)[counter] = var;
668  sumofstarts += starttime;
669  mindelta = MIN(mindelta, endtime - curtime); /* this amount of schifting holds for lb and ub */
670  counter++;
671  }
672
673  startindex--;
674  }
675
676  assert(mindelta > 0);
677  *lhs = lower ? sumofstarts + mindelta : sumofstarts - mindelta;
678
679  return SCIP_OKAY;
680 }
681
682 /** initialize the sorted event point arrays */
683 static
685  SCIP* scip, /**< SCIP data structure */
686  int nvars, /**< number of start time variables (activities) */
687  SCIP_VAR** vars, /**< array of start time variables */
688  int* durations, /**< array of durations per start time variable */
689  int* starttimes, /**< array to store sorted start events */
690  int* endtimes, /**< array to store sorted end events */
691  int* startindices, /**< permutation with rspect to the start times */
692  int* endindices, /**< permutation with rspect to the end times */
693  SCIP_Bool local /**< shall local bounds be used */
694  )
695 {
696  SCIP_VAR* var;
697  int j;
698
699  assert(vars != NULL || nvars == 0);
700
701  /* assign variables, start and endpoints to arrays */
702  for ( j = 0; j < nvars; ++j )
703  {
704  assert(vars != NULL);
705
706  var = vars[j];
707  assert(var != NULL);
708
709  if( local )
710  starttimes[j] = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
711  else
712  starttimes[j] = SCIPconvertRealToInt(scip, SCIPvarGetLbGlobal(var));
713
714  startindices[j] = j;
715
716  if( local )
717  endtimes[j] = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var)) + durations[j];
718  else
719  endtimes[j] = SCIPconvertRealToInt(scip, SCIPvarGetUbGlobal(var)) + durations[j];
720
721  endindices[j] = j;
722  }
723
724  /* sort the arrays not-decreasing according to startsolvalues and endsolvalues (and sort the indices in the same way) */
725  SCIPsortIntInt(starttimes, startindices, j);
726  SCIPsortIntInt(endtimes, endindices, j);
727 }
728
729 /** initialize the sorted event point arrays w.r.t. the given primal solutions */
730 static
732  SCIP* scip, /**< SCIP data structure */
733  SCIP_SOL* sol, /**< solution */
734  int nvars, /**< number of start time variables (activities) */
735  SCIP_VAR** vars, /**< array of start time variables */
736  int* durations, /**< array of durations per start time variable */
737  int* starttimes, /**< array to store sorted start events */
738  int* endtimes, /**< array to store sorted end events */
739  int* startindices, /**< permutation with rspect to the start times */
740  int* endindices /**< permutation with rspect to the end times */
741  )
742 {
743  SCIP_VAR* var;
744  int j;
745
746  assert(vars != NULL || nvars == 0);
747
748  /* assign variables, start and endpoints to arrays */
749  for ( j = 0; j < nvars; ++j )
750  {
751  assert(vars != NULL);
752
753  var = vars[j];
754  assert(var != NULL);
755
756  starttimes[j] = SCIPconvertRealToInt(scip, SCIPgetSolVal(scip, sol, var));
757  startindices[j] = j;
758
759  endtimes[j] = SCIPconvertRealToInt(scip, SCIPgetSolVal(scip, sol, var)) + durations[j];
760  endindices[j] = j;
761  }
762
763  /* sort the arrays not-decreasing according to startsolvalues and endsolvalues (and sort the indices in the same way) */
764  SCIPsortIntInt(starttimes, startindices, j);
765  SCIPsortIntInt(endtimes, endindices, j);
766 }
767
768 /** initialize the sorted event point arrays
769  *
770  * @todo Check the separation process!
771  */
772 static
774  SCIP* scip, /**< SCIP data structure */
775  SCIP_CONSDATA* consdata, /**< constraint data */
776  SCIP_SOL* sol, /**< primal CIP solution, NULL for current LP solution */
777  int* starttimes, /**< array to store sorted start events */
778  int* endtimes, /**< array to store sorted end events */
779  int* startindices, /**< permutation with rspect to the start times */
780  int* endindices, /**< permutation with rspect to the end times */
781  int* nvars, /**< number of variables that are integral */
782  SCIP_Bool lower /**< shall the constraints be derived for lower or upper bounds? */
783  )
784 {
785  SCIP_VAR* var;
786  int tmpnvars;
787  int j;
788
789  tmpnvars = consdata->nvars;
790  *nvars = 0;
791
792  /* assign variables, start and endpoints to arrays */
793  for ( j = 0; j < tmpnvars; ++j )
794  {
795  var = consdata->vars[j];
796  assert(var != NULL);
797  assert(consdata->durations[j] > 0);
798  assert(consdata->demands[j] > 0);
799
800  if( lower )
801  {
802  /* only consider jobs that are at their lower or upper bound */
803  if( !SCIPisFeasIntegral(scip, SCIPgetSolVal(scip, sol, var))
804  || !SCIPisFeasEQ(scip, SCIPgetSolVal(scip, sol, var), SCIPvarGetLbLocal(var)) )
805  continue;
806
807  starttimes[*nvars] = SCIPconvertRealToInt(scip, SCIPgetSolVal(scip, sol, var));
808  startindices[*nvars] = j;
809
810  endtimes[*nvars] = starttimes[*nvars] + consdata->durations[j];
811  endindices[*nvars] = j;
812
813  SCIPdebugMsg(scip, "%d: variable <%s>[%g,%g] (sol %g, duration %d) starttime %d, endtime = %d, demand = %d\n",
814  *nvars, SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), SCIPgetSolVal(scip, sol, var),
815  consdata->durations[j],
816  starttimes[*nvars], starttimes[*nvars] + consdata->durations[startindices[*nvars]],
817  consdata->demands[startindices[*nvars]]);
818
819  (*nvars)++;
820  }
821  else
822  {
823  if( !SCIPisFeasIntegral(scip, SCIPgetSolVal(scip, sol, var))
824  || !SCIPisFeasEQ(scip, SCIPgetSolVal(scip, sol, var), SCIPvarGetUbLocal(var)) )
825  continue;
826
827  starttimes[*nvars] = SCIPconvertRealToInt(scip, SCIPgetSolVal(scip, sol, var));
828  startindices[*nvars] = j;
829
830  endtimes[*nvars] = starttimes[*nvars] + consdata->durations[j];
831  endindices[*nvars] = j;
832
833  SCIPdebugMsg(scip, "%d: variable <%s>[%g,%g] (sol %g, duration %d) starttime %d, endtime = %d, demand = %d\n",
834  *nvars, SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), SCIPgetSolVal(scip, sol, var),
835  consdata->durations[j],
836  starttimes[*nvars], starttimes[*nvars] + consdata->durations[startindices[*nvars]],
837  consdata->demands[startindices[*nvars]]);
838
839  (*nvars)++;
840  }
841  }
842
843  /* sort the arrays not-decreasing according to startsolvalues and endsolvalues (and sort the indices in the same way) */
844  SCIPsortIntInt(starttimes, startindices, *nvars);
845  SCIPsortIntInt(endtimes, endindices, *nvars);
846
847 #ifdef SCIP_DEBUG
848  SCIPdebugMsg(scip, "sorted output %d\n", *nvars);
849
850  for ( j = 0; j < *nvars; ++j )
851  {
852  SCIPdebugMsg(scip, "%d: job[%d] starttime %d, endtime = %d, demand = %d\n", j,
853  startindices[j], starttimes[j], starttimes[j] + consdata->durations[startindices[j]],
854  consdata->demands[startindices[j]]);
855  }
856
857  for ( j = 0; j < *nvars; ++j )
858  {
859  SCIPdebugMsg(scip, "%d: job[%d] endtime %d, demand = %d\n", j, endindices[j], endtimes[j],
860  consdata->demands[endindices[j]]);
861  }
862 #endif
863 }
864
865 #ifdef SCIP_STATISTIC
866 /** this method checks for relevant intervals for energetic reasoning */
867 static
868 SCIP_RETCODE computeRelevantEnergyIntervals(
869  SCIP* scip, /**< SCIP data structure */
870  int nvars, /**< number of start time variables (activities) */
871  SCIP_VAR** vars, /**< array of start time variables */
872  int* durations, /**< array of durations */
873  int* demands, /**< array of demands */
874  int capacity, /**< cumulative capacity */
875  int hmin, /**< left bound of time axis to be considered (including hmin) */
876  int hmax, /**< right bound of time axis to be considered (not including hmax) */
877  int** timepoints, /**< array to store relevant points in time */
878  SCIP_Real** cumulativedemands, /**< array to store the estimated cumulative demand for each point in time */
879  int* ntimepoints, /**< pointer to store the number of timepoints */
880  int* maxdemand, /**< pointer to store maximum over all demands */
881  SCIP_Real* minfreecapacity /**< pointer to store the minimum free capacity */
882  )
883 {
884  int* starttimes; /* stores when each job is starting */
885  int* endtimes; /* stores when each job ends */
886  int* startindices; /* we will sort the startsolvalues, thus we need to know wich index of a job it corresponds to */
887  int* endindices; /* we will sort the endsolvalues, thus we need to know wich index of a job it corresponds to */
888
889  SCIP_Real totaldemand;
890  int curtime; /* point in time which we are just checking */
891  int endindex; /* index of endsolvalues with: endsolvalues[endindex] > curtime */
892
893  int j;
894
895  assert( scip != NULL );
896  assert(durations != NULL);
897  assert(demands != NULL);
898  assert(capacity >= 0);
899
900  /* if no activities are associated with this cumulative then this constraint is redundant */
901  if( nvars == 0 )
902  return SCIP_OKAY;
903
904  assert(vars != NULL);
905
906  SCIP_CALL( SCIPallocBufferArray(scip, &starttimes, nvars) );
907  SCIP_CALL( SCIPallocBufferArray(scip, &endtimes, nvars) );
908  SCIP_CALL( SCIPallocBufferArray(scip, &startindices, nvars) );
909  SCIP_CALL( SCIPallocBufferArray(scip, &endindices, nvars) );
910
911  /* create event point arrays */
912  createSortedEventpoints(scip, nvars, vars, durations, starttimes, endtimes, startindices, endindices, TRUE);
913
914  endindex = 0;
915  totaldemand = 0.0;
916
917  *ntimepoints = 0;
918  (*timepoints)[0] = starttimes[0];
919  (*cumulativedemands)[0] = 0;
920  *maxdemand = 0;
921
922  /* check each startpoint of a job whether the capacity is kept or not */
923  for( j = 0; j < nvars; ++j )
924  {
925  int lct;
926  int idx;
927
928  curtime = starttimes[j];
929
930  if( curtime >= hmax )
931  break;
932
933  /* free all capacity usages of jobs the are no longer running */
934  while( endindex < nvars && endtimes[endindex] <= curtime )
935  {
936  int est;
937
938  if( (*timepoints)[*ntimepoints] < endtimes[endindex] )
939  {
940  (*ntimepoints)++;
941  (*timepoints)[*ntimepoints] = endtimes[endindex];
942  (*cumulativedemands)[*ntimepoints] = 0;
943  }
944
945  idx = endindices[endindex];
946  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[idx]));
947  totaldemand -= (SCIP_Real) demands[idx] * durations[idx] / (endtimes[endindex] - est);
948  endindex++;
949
950  (*cumulativedemands)[*ntimepoints] = totaldemand;
951  }
952
953  idx = startindices[j];
954  lct = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[idx]) + durations[idx]);
955  totaldemand += (SCIP_Real) demands[idx] * durations[idx] / (lct - starttimes[j]);
956
957  if( (*timepoints)[*ntimepoints] < curtime )
958  {
959  (*ntimepoints)++;
960  (*timepoints)[*ntimepoints] = curtime;
961  (*cumulativedemands)[*ntimepoints] = 0;
962  }
963
964  (*cumulativedemands)[*ntimepoints] = totaldemand;
965
966  /* add the relative capacity requirements for all job which start at the curtime */
967  while( j+1 < nvars && starttimes[j+1] == curtime )
968  {
969  ++j;
970  idx = startindices[j];
971  lct = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[idx]) + durations[idx]);
972  totaldemand += (SCIP_Real) demands[idx] * durations[idx] / (lct - starttimes[j]);
973
974  (*cumulativedemands)[*ntimepoints] = totaldemand;
975  }
976  } /*lint --e{850}*/
977
978  /* free all capacity usages of jobs that are no longer running */
979  while( endindex < nvars/* && endtimes[endindex] < hmax*/)
980  {
981  int est;
982  int idx;
983
984  if( (*timepoints)[*ntimepoints] < endtimes[endindex] )
985  {
986  (*ntimepoints)++;
987  (*timepoints)[*ntimepoints] = endtimes[endindex];
988  (*cumulativedemands)[*ntimepoints] = 0;
989  }
990
991  idx = endindices[endindex];
992  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[idx]));
993  totaldemand -= (SCIP_Real) demands[idx] * durations[idx] / (endtimes[endindex] - est);
994  (*cumulativedemands)[*ntimepoints] = totaldemand;
995
996  ++endindex;
997  }
998
999  (*ntimepoints)++;
1000  /* compute minimum free capacity */
1001  (*minfreecapacity) = INT_MAX;
1002  for( j = 0; j < *ntimepoints; ++j )
1003  {
1004  if( (*timepoints)[j] >= hmin && (*timepoints)[j] < hmax )
1005  *minfreecapacity = MIN( *minfreecapacity, (SCIP_Real)capacity - (*cumulativedemands)[j] );
1006  }
1007
1008  /* free buffer arrays */
1009  SCIPfreeBufferArray(scip, &endindices);
1010  SCIPfreeBufferArray(scip, &startindices);
1011  SCIPfreeBufferArray(scip, &endtimes);
1012  SCIPfreeBufferArray(scip, &starttimes);
1013
1014  return SCIP_OKAY;
1015 }
1016
1017 /** evaluates the cumulativeness and disjointness factor of a cumulative constraint */
1018 static
1019 SCIP_RETCODE evaluateCumulativeness(
1020  SCIP* scip, /**< pointer to scip */
1021  SCIP_CONS* cons /**< cumulative constraint */
1022  )
1023 {
1024  SCIP_CONSDATA* consdata;
1025  int nvars;
1026  int v;
1027  int capacity;
1028
1029  /* output values: */
1030  SCIP_Real disjfactor2; /* (peak-capacity)/capacity * (large demands/nvars_t) */
1031  SCIP_Real cumfactor1;
1032  SCIP_Real resstrength1; /* overall strength */
1033  SCIP_Real resstrength2; /* timepoint wise maximum */
1034
1035  /* helpful variables: */
1036  SCIP_Real globalpeak;
1037  SCIP_Real globalmaxdemand;
1038
1039  /* get constraint data structure */
1040  consdata = SCIPconsGetData(cons);
1041  assert(consdata != NULL);
1042
1043  nvars = consdata->nvars;
1044  capacity = consdata->capacity;
1045  globalpeak = 0.0;
1046  globalmaxdemand = 0.0;
1047
1048  disjfactor2 = 0.0;
1049  cumfactor1 = 0.0;
1050  resstrength2 = 0.0;
1051
1052  /* check each starting time (==each job, but inefficient) */
1053  for( v = 0; v < nvars; ++v )
1054  {
1055  SCIP_Real peak;
1056  SCIP_Real maxdemand;
1058  int ndemands;
1059  int nlarge;
1060
1061  int timepoint;
1062  int j;
1063  timepoint = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(consdata->vars[v]));
1064  peak = consdata->demands[v];
1065  ndemands = 1;
1066  maxdemand = 0;
1067  nlarge = 0;
1068
1069  if( consdata->demands[v] > capacity / 3 )
1070  nlarge++;
1071
1072  for( j = 0; j < nvars; ++j )
1073  {
1074  int lb;
1075
1076  if( j == v )
1077  continue;
1078
1079  maxdemand = 0.0;
1080  lb = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(consdata->vars[j]));
1081
1082  if( lb <= timepoint && lb + consdata->durations[j] > timepoint )
1083  {
1084  peak += consdata->demands[j];
1085  ndemands++;
1086
1087  if( consdata->demands[j] > consdata->capacity / 3 )
1088  nlarge++;
1089  }
1090  }
1091
1092  deltademand = (SCIP_Real)peak / (SCIP_Real)ndemands;
1093  globalpeak = MAX(globalpeak, peak);
1094  globalmaxdemand = MAX(globalmaxdemand, maxdemand);
1095
1096  if( peak > capacity )
1097  {
1098  disjfactor2 = MAX( disjfactor2, (peak-(SCIP_Real)capacity)/peak * (nlarge/(SCIP_Real)ndemands) );
1099  cumfactor1 = MAX( cumfactor1, (peak-capacity)/peak * (capacity-deltademand)/(SCIP_Real)capacity );
1100  resstrength2 = MAX(resstrength2, (capacity-maxdemand)/(peak-maxdemand) );
1101  }
1102  }
1103
1104  resstrength1 = (capacity-globalmaxdemand) / (globalpeak-globalmaxdemand);
1105
1106  consdata->maxpeak = SCIPconvertRealToInt(scip, globalpeak);
1107  consdata->disjfactor2 = disjfactor2;
1108  consdata->cumfactor1 = cumfactor1;
1109  consdata->resstrength2 = resstrength2;
1110  consdata->resstrength1 = resstrength1;
1111
1112  /* get estimated res strength */
1113  {
1114  int* timepoints;
1115  SCIP_Real* estimateddemands;
1116  int ntimepoints;
1117  int maxdemand;
1118  SCIP_Real minfreecapacity;
1119
1120  SCIP_CALL( SCIPallocBufferArray(scip, &timepoints, 2*nvars) );
1121  SCIP_CALL( SCIPallocBufferArray(scip, &estimateddemands, 2*nvars) );
1122
1123  ntimepoints = 0;
1124  minfreecapacity = INT_MAX;
1125
1126  SCIP_CALL( computeRelevantEnergyIntervals(scip, nvars, consdata->vars,
1127  consdata->durations, consdata->demands,
1128  capacity, consdata->hmin, consdata->hmax, &timepoints, &estimateddemands,
1129  &ntimepoints, &maxdemand, &minfreecapacity) );
1130
1131  /* free buffer arrays */
1132  SCIPfreeBufferArray(scip, &estimateddemands);
1133  SCIPfreeBufferArray(scip, &timepoints);
1134
1135  consdata->estimatedstrength = (SCIP_Real)(capacity - minfreecapacity) / (SCIP_Real) capacity;
1136  }
1137
1138  SCIPstatisticPrintf("cumulative constraint<%s>: DISJ1=%g, DISJ2=%g, CUM=%g, RS1 = %g, RS2 = %g, EST = %g\n",
1139  SCIPconsGetName(cons), consdata->disjfactor1, disjfactor2, cumfactor1, resstrength1, resstrength2,
1140  consdata->estimatedstrength);
1141
1142  return SCIP_OKAY;
1143 }
1144 #endif
1145
1146 /** gets the active variables together with the constant */
1147 static
1149  SCIP* scip, /**< SCIP data structure */
1150  SCIP_VAR** var, /**< pointer to store the active variable */
1151  int* scalar, /**< pointer to store the scalar */
1152  int* constant /**< pointer to store the constant */
1153  )
1154 {
1155  if( !SCIPvarIsActive(*var) )
1156  {
1157  SCIP_Real realscalar;
1158  SCIP_Real realconstant;
1159
1160  realscalar = 1.0;
1161  realconstant = 0.0;
1162
1164
1165  /* transform variable to active variable */
1166  SCIP_CALL( SCIPgetProbvarSum(scip, var, &realscalar, &realconstant) );
1167  assert(!SCIPisZero(scip, realscalar));
1168  assert(SCIPvarIsActive(*var));
1169
1170  if( realconstant < 0.0 )
1171  (*constant) = -SCIPconvertRealToInt(scip, -realconstant);
1172  else
1173  (*constant) = SCIPconvertRealToInt(scip, realconstant);
1174
1175  if( realscalar < 0.0 )
1176  (*scalar) = -SCIPconvertRealToInt(scip, -realscalar);
1177  else
1178  (*scalar) = SCIPconvertRealToInt(scip, realscalar);
1179  }
1180  else
1181  {
1182  (*scalar) = 1;
1183  (*constant) = 0;
1184  }
1185
1186  assert(*scalar != 0);
1187
1188  return SCIP_OKAY;
1189 }
1190
1191 /** computes the total energy of all jobs */
1192 static
1194  int* durations, /**< array of job durations */
1195  int* demands, /**< array of job demands */
1196  int njobs /**< number of jobs */
1197  )
1198 {
1199  SCIP_Longint energy;
1200  int j;
1201
1202  energy = 0;
1203
1204  for( j = 0; j < njobs; ++j )
1205  energy += (SCIP_Longint) durations[j] * demands[j];
1206
1207  return energy;
1208 }
1209
1210 /**@} */
1211
1212 /**@name Default method to solve a cumulative condition
1213  *
1214  * @{
1215  */
1216
1217 /** setup and solve subscip to solve single cumulative condition */
1218 static
1220  SCIP* subscip, /**< subscip data structure */
1221  SCIP_Real* objvals, /**< array of objective coefficients for each job (linear objective function), or NULL if none */
1222  int* durations, /**< array of durations */
1223  int* demands, /**< array of demands */
1224  int njobs, /**< number of jobs (activities) */
1225  int capacity, /**< cumulative capacity */
1226  int hmin, /**< left bound of time axis to be considered (including hmin) */
1227  int hmax, /**< right bound of time axis to be considered (not including hmax) */
1228  SCIP_Longint maxnodes, /**< maximum number of branch-and-bound nodes (-1: no limit) */
1229  SCIP_Real timelimit, /**< time limit for solving in seconds */
1230  SCIP_Real memorylimit, /**< memory limit for solving in mega bytes (MB) */
1231  SCIP_Real* ests, /**< array of earliest start times for each job */
1232  SCIP_Real* lsts, /**< array of latest start times for each job */
1233  SCIP_Bool* infeasible, /**< pointer to store if the subproblem was infeasible */
1234  SCIP_Bool* unbounded, /**< pointer to store if the problem is unbounded */
1235  SCIP_Bool* solved, /**< pointer to store if the problem is solved (to optimality) */
1236  SCIP_Bool* error /**< pointer to store if an error occurred */
1237  )
1238 {
1239  SCIP_VAR** subvars;
1240  SCIP_CONS* cons;
1241
1242  char name[SCIP_MAXSTRLEN];
1243  int v;
1244  SCIP_RETCODE retcode;
1245
1246  assert(subscip != NULL);
1247
1248  /* copy all plugins */
1250
1251  /* create the subproblem */
1252  SCIP_CALL( SCIPcreateProbBasic(subscip, "cumulative") );
1253
1254  SCIP_CALL( SCIPallocBlockMemoryArray(subscip, &subvars, njobs) );
1255
1256  /* create for each job a start time variable */
1257  for( v = 0; v < njobs; ++v )
1258  {
1259  SCIP_Real objval;
1260
1261  /* construct variable name */
1262  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "job%d", v);
1263
1264  if( objvals == NULL )
1265  objval = 0.0;
1266  else
1267  objval = objvals[v];
1268
1269  SCIP_CALL( SCIPcreateVarBasic(subscip, &subvars[v], name, ests[v], lsts[v], objval, SCIP_VARTYPE_INTEGER) );
1270  SCIP_CALL( SCIPaddVar(subscip, subvars[v]) );
1271  }
1272
1273  /* create cumulative constraint */
1274  SCIP_CALL( SCIPcreateConsBasicCumulative(subscip, &cons, "cumulative",
1275  njobs, subvars, durations, demands, capacity) );
1276
1277  /* set effective horizon */
1278  SCIP_CALL( SCIPsetHminCumulative(subscip, cons, hmin) );
1279  SCIP_CALL( SCIPsetHmaxCumulative(subscip, cons, hmax) );
1280
1281  /* add cumulative constraint */
1282  SCIP_CALL( SCIPaddCons(subscip, cons) );
1283  SCIP_CALL( SCIPreleaseCons(subscip, &cons) );
1284
1285  /* set CP solver settings
1286  *
1287  * @note This "meta" setting has to be set first since this call overwrite all parameters including for example the
1288  * time limit.
1289  */
1291
1292  /* do not abort subproblem on CTRL-C */
1293  SCIP_CALL( SCIPsetBoolParam(subscip, "misc/catchctrlc", FALSE) );
1294
1295  /* disable output to console */
1296  SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
1297
1298  /* set limits for the subproblem */
1299  SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", maxnodes) );
1300  SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", timelimit) );
1301  SCIP_CALL( SCIPsetRealParam(subscip, "limits/memory", memorylimit) );
1302
1303  /* forbid recursive call of heuristics and separators solving subMIPs */
1304  SCIP_CALL( SCIPsetSubscipsOff(subscip, TRUE) );
1305
1306  /* solve single cumulative constraint by branch and bound */
1307  retcode = SCIPsolve(subscip);
1308
1309  if( retcode != SCIP_OKAY )
1310  (*error) = TRUE;
1311  else
1312  {
1313  SCIPdebugMsg(subscip, "solved single cumulative condition with status %d\n", SCIPgetStatus(subscip));
1314
1315  /* evaluated solution status */
1316  switch( SCIPgetStatus(subscip) )
1317  {
1318  case SCIP_STATUS_INFORUNBD:
1320  (*infeasible) = TRUE;
1321  (*solved) = TRUE;
1322  break;
1323  case SCIP_STATUS_UNBOUNDED:
1324  (*unbounded) = TRUE;
1325  (*solved) = TRUE;
1326  break;
1327  case SCIP_STATUS_OPTIMAL:
1328  {
1329  SCIP_SOL* sol;
1330  SCIP_Real solval;
1331
1332  sol = SCIPgetBestSol(subscip);
1333  assert(sol != NULL);
1334
1335  for( v = 0; v < njobs; ++v )
1336  {
1337  solval = SCIPgetSolVal(subscip, sol, subvars[v]);
1338
1339  ests[v] = solval;
1340  lsts[v] = solval;
1341  }
1342  (*solved) = TRUE;
1343  break;
1344  }
1345  case SCIP_STATUS_NODELIMIT:
1347  case SCIP_STATUS_TIMELIMIT:
1348  case SCIP_STATUS_MEMLIMIT:
1350  case SCIP_STATUS_TERMINATE:
1351  /* transfer the global bound changes */
1352  for( v = 0; v < njobs; ++v )
1353  {
1354  ests[v] = SCIPvarGetLbGlobal(subvars[v]);
1355  lsts[v] = SCIPvarGetUbGlobal(subvars[v]);
1356  }
1357  (*solved) = FALSE;
1358  break;
1359
1360  case SCIP_STATUS_UNKNOWN:
1362  case SCIP_STATUS_GAPLIMIT:
1363  case SCIP_STATUS_SOLLIMIT:
1366  SCIPerrorMessage("invalid status code <%d>\n", SCIPgetStatus(subscip));
1367  return SCIP_INVALIDDATA;
1368  }
1369  }
1370
1371  /* release all variables */
1372  for( v = 0; v < njobs; ++v )
1373  {
1374  SCIP_CALL( SCIPreleaseVar(subscip, &subvars[v]) );
1375  }
1376
1377  SCIPfreeBlockMemoryArray(subscip, &subvars, njobs);
1378
1379  return SCIP_OKAY;
1380 }
1381
1382 /** solve single cumulative condition using SCIP and a single cumulative constraint */
1383 static
1384 SCIP_DECL_SOLVECUMULATIVE(solveCumulativeViaScipCp)
1386  SCIP* subscip;
1387
1388  SCIP_RETCODE retcode;
1389
1390  assert(njobs > 0);
1391
1392  (*solved) = FALSE;
1393  (*infeasible) = FALSE;
1394  (*unbounded) = FALSE;
1395  (*error) = FALSE;
1396
1397  SCIPdebugMessage("solve independent cumulative condition with %d variables\n", njobs);
1398
1399  /* initialize the sub-problem */
1400  SCIP_CALL( SCIPcreate(&subscip) );
1401
1402  /* create and solve the subproblem. catch possible errors */
1403  retcode = setupAndSolveCumulativeSubscip(subscip, objvals, durations, demands,
1404  njobs, capacity, hmin, hmax,
1405  maxnodes, timelimit, memorylimit,
1406  ests, lsts,
1407  infeasible, unbounded, solved, error);
1408
1409  /* free the subscip in any case */
1410  SCIP_CALL( SCIPfree(&subscip) );
1411
1412  SCIP_CALL( retcode );
1413
1414  return SCIP_OKAY;
1415 }
1416
1417 #if 0
1418 /** solve single cumulative condition using SCIP and the time indexed formulation */
1419 static
1420 SCIP_DECL_SOLVECUMULATIVE(solveCumulativeViaScipMip)
1421 {
1422  SCIP* subscip;
1423  SCIP_VAR*** binvars;
1424  SCIP_RETCODE retcode;
1425  char name[SCIP_MAXSTRLEN];
1426  int minest;
1427  int maxlct;
1428  int t;
1429  int v;
1430
1431  assert(njobs > 0);
1432
1433  (*solved) = FALSE;
1434  (*infeasible) = FALSE;
1435  (*unbounded) = FALSE;
1436  (*error) = FALSE;
1437
1438  SCIPdebugMsg(scip, "solve independent cumulative condition with %d variables\n", njobs);
1439
1440  /* initialize the sub-problem */
1441  SCIP_CALL( SCIPcreate(&subscip) );
1442
1443  /* copy all plugins */
1445
1446  /* create the subproblem */
1447  SCIP_CALL( SCIPcreateProbBasic(subscip, "cumulative") );
1448
1449  SCIP_CALL( SCIPallocBufferArray(subscip, &binvars, njobs) );
1450
1451  minest = INT_MAX;
1452  maxlct = INT_MIN;
1453
1454  /* create for each job and time step a binary variable which is one if this jobs starts at this time point and a set
1455  * partitioning constrain which forces that job starts
1456  */
1457  for( v = 0; v < njobs; ++v )
1458  {
1459  SCIP_CONS* cons;
1460  SCIP_Real objval;
1461  int timeinterval;
1462  int est;
1463  int lst;
1464
1465  if( objvals == NULL )
1466  objval = 0.0;
1467  else
1468  objval = objvals[v];
1469
1470  est = ests[v];
1471  lst = lsts[v];
1472
1473  /* compute number of possible start points */
1474  timeinterval = lst - est + 1;
1475  assert(timeinterval > 0);
1476
1477  /* compute the smallest earliest start time and largest latest completion time */
1478  minest = MIN(minest, est);
1479  maxlct = MAX(maxlct, lst + durations[v]);
1480
1481  /* construct constraint name */
1482  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "job_%d", v);
1483
1484  SCIP_CALL( SCIPcreateConsBasicSetpart(subscip, &cons, name, 0, NULL) );
1485
1486  SCIP_CALL( SCIPallocBufferArray(subscip, &binvars[v], timeinterval) );
1487
1488  for( t = 0; t < timeinterval; ++t )
1489  {
1490  SCIP_VAR* binvar;
1491
1492  /* construct varibale name */
1493  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "job_%d_time_%d", v, t + est);
1494
1495  SCIP_CALL( SCIPcreateVarBasic(subscip, &binvar, name, 0.0, 1.0, objval, SCIP_VARTYPE_BINARY) );
1496  SCIP_CALL( SCIPaddVar(subscip, binvar) );
1497
1498  /* add binary varibale to the set partitioning constraint which ensures that the job is started */
1499  SCIP_CALL( SCIPaddCoefSetppc(subscip, cons, binvar) );
1500
1501  binvars[v][t] = binvar;
1502  }
1503
1504  /* add and release the set partitioning constraint */
1505  SCIP_CALL( SCIPaddCons(subscip, cons) );
1506  SCIP_CALL( SCIPreleaseCons(subscip, &cons) );
1507  }
1508
1509  /* adjusted the smallest earliest start time and the largest latest completion time with the effective horizon */
1510  hmin = MAX(hmin, minest);
1511  hmax = MIN(hmax, maxlct);
1512  assert(hmin > INT_MIN);
1513  assert(hmax < INT_MAX);
1514  assert(hmin < hmax);
1515
1516  /* create for each time a knapsack constraint which ensures that the resource capacity is not exceeded */
1517  for( t = hmin; t < hmax; ++t )
1518  {
1519  SCIP_CONS* cons;
1520
1521  /* construct constraint name */
1522  (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "time_%d", t);
1523
1524  /* create an empty knapsack constraint */
1525  SCIP_CALL( SCIPcreateConsBasicKnapsack(subscip, &cons, name, 0, NULL, NULL, (SCIP_Longint)capacity) );
1526
1527  /* add all jobs which potentially can be processed at that time point */
1528  for( v = 0; v < njobs; ++v )
1529  {
1530  int duration;
1531  int demand;
1532  int start;
1533  int end;
1534  int est;
1535  int lst;
1536  int k;
1537
1538  est = ests[v];
1539  lst = lsts[v] ;
1540
1541  duration = durations[v];
1542  assert(duration > 0);
1543
1544  /* check if the varibale is processed potentially at time point t */
1545  if( t < est || t >= lst + duration )
1546  continue;
1547
1548  demand = demands[v];
1549  assert(demand >= 0);
1550
1551  start = MAX(t - duration + 1, est);
1552  end = MIN(t, lst);
1553
1554  assert(start <= end);
1555
1556  for( k = start; k <= end; ++k )
1557  {
1558  assert(binvars[v][k] != NULL);
1559  SCIP_CALL( SCIPaddCoefKnapsack(subscip, cons, binvars[v][k], (SCIP_Longint) demand) );
1560  }
1561  }
1562
1563  /* add and release the knapsack constraint */
1564  SCIP_CALL( SCIPaddCons(subscip, cons) );
1565  SCIP_CALL( SCIPreleaseCons(subscip, &cons) );
1566  }
1567
1568  /* do not abort subproblem on CTRL-C */
1569  SCIP_CALL( SCIPsetBoolParam(subscip, "misc/catchctrlc", FALSE) );
1570
1571  /* disable output to console */
1572  SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
1573
1574  /* set limits for the subproblem */
1575  SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", maxnodes) );
1576  SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", timelimit) );
1577  SCIP_CALL( SCIPsetRealParam(subscip, "limits/memory", memorylimit) );
1578
1579  /* solve single cumulative constraint by branch and bound */
1580  retcode = SCIPsolve(subscip);
1581
1582  if( retcode != SCIP_OKAY )
1583  (*error) = TRUE;
1584  else
1585  {
1586  SCIPdebugMsg(scip, "solved single cumulative condition with status %d\n", SCIPgetStatus(subscip));
1587
1588  /* evaluated solution status */
1589  switch( SCIPgetStatus(subscip) )
1590  {
1591  case SCIP_STATUS_INFORUNBD:
1593  (*infeasible) = TRUE;
1594  (*solved) = TRUE;
1595  break;
1596  case SCIP_STATUS_UNBOUNDED:
1597  (*unbounded) = TRUE;
1598  (*solved) = TRUE;
1599  break;
1600  case SCIP_STATUS_OPTIMAL:
1601  {
1602  SCIP_SOL* sol;
1603
1604  sol = SCIPgetBestSol(subscip);
1605  assert(sol != NULL);
1606
1607  for( v = 0; v < njobs; ++v )
1608  {
1609  int timeinterval;
1610  int est;
1611  int lst;
1612
1613  est = ests[v];
1614  lst = lsts[v];
1615
1616  /* compute number of possible start points */
1617  timeinterval = lst - est + 1;
1618
1619  /* check which binary varibale is set to one */
1620  for( t = 0; t < timeinterval; ++t )
1621  {
1622  if( SCIPgetSolVal(subscip, sol, binvars[v][t]) > 0.5 )
1623  {
1624  ests[v] = est + t;
1625  lsts[v] = est + t;
1626  break;
1627  }
1628  }
1629  }
1630
1631  (*solved) = TRUE;
1632  break;
1633  }
1634  case SCIP_STATUS_NODELIMIT:
1636  case SCIP_STATUS_TIMELIMIT:
1637  case SCIP_STATUS_MEMLIMIT:
1639  /* transfer the global bound changes */
1640  for( v = 0; v < njobs; ++v )
1641  {
1642  int timeinterval;
1643  int est;
1644  int lst;
1645
1646  est = ests[v];
1647  lst = lsts[v];
1648
1649  /* compute number of possible start points */
1650  timeinterval = lst - est + 1;
1651
1652  /* check which binary varibale is the first binary varibale which is not globally fixed to zero */
1653  for( t = 0; t < timeinterval; ++t )
1654  {
1655  if( SCIPvarGetUbGlobal(binvars[v][t]) > 0.5 )
1656  {
1657  ests[v] = est + t;
1658  break;
1659  }
1660  }
1661
1662  /* check which binary varibale is the last binary varibale which is not globally fixed to zero */
1663  for( t = timeinterval - 1; t >= 0; --t )
1664  {
1665  if( SCIPvarGetUbGlobal(binvars[v][t]) > 0.5 )
1666  {
1667  lsts[v] = est + t;
1668  break;
1669  }
1670  }
1671  }
1672  (*solved) = FALSE;
1673  break;
1674
1675  case SCIP_STATUS_UNKNOWN:
1677  case SCIP_STATUS_GAPLIMIT:
1678  case SCIP_STATUS_SOLLIMIT:
1680  SCIPerrorMessage("invalid status code <%d>\n", SCIPgetStatus(subscip));
1681  return SCIP_INVALIDDATA;
1682  }
1683  }
1684
1685  /* release all variables */
1686  for( v = 0; v < njobs; ++v )
1687  {
1688  int timeinterval;
1689  int est;
1690  int lst;
1691
1692  est = ests[v];
1693  lst = lsts[v];
1694
1695  /* compute number of possible start points */
1696  timeinterval = lst - est + 1;
1697
1698  for( t = 0; t < timeinterval; ++t )
1699  {
1700  SCIP_CALL( SCIPreleaseVar(subscip, &binvars[v][t]) );
1701  }
1702  SCIPfreeBufferArray(subscip, &binvars[v]);
1703  }
1704
1705  SCIPfreeBufferArray(subscip, &binvars);
1706
1707  SCIP_CALL( SCIPfree(&subscip) );
1708
1709  return SCIP_OKAY;
1710 }
1711 #endif
1712
1713 /**@} */
1714
1715 /**@name Constraint handler data
1716  *
1717  * Method used to create and free the constraint handler data when including and removing the cumulative constraint
1718  * handler.
1719  *
1720  * @{
1721  */
1722
1723 /** creates constaint handler data for cumulative constraint handler */
1724 static
1726  SCIP* scip, /**< SCIP data structure */
1727  SCIP_CONSHDLRDATA** conshdlrdata, /**< pointer to store the constraint handler data */
1728  SCIP_EVENTHDLR* eventhdlr /**< event handler */
1729  )
1730 {
1731  /* create precedence constraint handler data */
1732  assert(scip != NULL);
1733  assert(conshdlrdata != NULL);
1734  assert(eventhdlr != NULL);
1735
1736  SCIP_CALL( SCIPallocBlockMemory(scip, conshdlrdata) );
1737
1738  /* set event handler for checking if bounds of start time variables are tighten */
1739  (*conshdlrdata)->eventhdlr = eventhdlr;
1740
1741  /* set default methed for solving single cumulative conditions using SCIP and a CP model */
1742  (*conshdlrdata)->solveCumulative = solveCumulativeViaScipCp;
1743
1744 #ifdef SCIP_STATISTIC
1745  (*conshdlrdata)->nlbtimetable = 0;
1746  (*conshdlrdata)->nubtimetable = 0;
1747  (*conshdlrdata)->ncutofftimetable = 0;
1748  (*conshdlrdata)->nlbedgefinder = 0;
1749  (*conshdlrdata)->nubedgefinder = 0;
1750  (*conshdlrdata)->ncutoffedgefinder = 0;
1751  (*conshdlrdata)->ncutoffoverload = 0;
1752  (*conshdlrdata)->ncutoffoverloadTTEF = 0;
1753
1754  (*conshdlrdata)->nirrelevantjobs = 0;
1755  (*conshdlrdata)->nalwaysruns = 0;
1756  (*conshdlrdata)->nremovedlocks = 0;
1757  (*conshdlrdata)->ndualfixs = 0;
1758  (*conshdlrdata)->ndecomps = 0;
1759  (*conshdlrdata)->ndualbranchs = 0;
1760  (*conshdlrdata)->nallconsdualfixs = 0;
1761  (*conshdlrdata)->naddedvarbounds = 0;
1762  (*conshdlrdata)->naddeddisjunctives = 0;
1763 #endif
1764
1765  return SCIP_OKAY;
1766 }
1767
1768 /** frees constraint handler data for logic or constraint handler */
1769 static
1770 void conshdlrdataFree(
1771  SCIP* scip, /**< SCIP data structure */
1772  SCIP_CONSHDLRDATA** conshdlrdata /**< pointer to the constraint handler data */
1773  )
1774 {
1775  assert(conshdlrdata != NULL);
1776  assert(*conshdlrdata != NULL);
1777
1778  SCIPfreeBlockMemory(scip, conshdlrdata);
1779 }
1780
1781 /**@} */
1782
1783
1784 /**@name Constraint data methods
1785  *
1786  * @{
1787  */
1788
1789 /** catches bound change events for all variables in transformed cumulative constraint */
1790 static
1792  SCIP* scip, /**< SCIP data structure */
1793  SCIP_CONSDATA* consdata, /**< cumulative constraint data */
1794  SCIP_EVENTHDLR* eventhdlr /**< event handler to call for the event processing */
1795  )
1796 {
1797  int v;
1798
1799  assert(scip != NULL);
1800  assert(consdata != NULL);
1801  assert(eventhdlr != NULL);
1802
1803  /* catch event for every single variable */
1804  for( v = 0; v < consdata->nvars; ++v )
1805  {
1806  SCIP_CALL( SCIPcatchVarEvent(scip, consdata->vars[v],
1807  SCIP_EVENTTYPE_BOUNDTIGHTENED, eventhdlr, (SCIP_EVENTDATA*)consdata, NULL) );
1808  }
1809
1810  return SCIP_OKAY;
1811 }
1812
1813 /** drops events for variable at given position */
1814 static
1816  SCIP* scip, /**< SCIP data structure */
1817  SCIP_CONSDATA* consdata, /**< cumulative constraint data */
1818  SCIP_EVENTHDLR* eventhdlr, /**< event handler to call for the event processing */
1819  int pos /**< array position of variable to catch bound change events for */
1820  )
1821 {
1822  assert(scip != NULL);
1823  assert(consdata != NULL);
1824  assert(eventhdlr != NULL);
1825  assert(0 <= pos && pos < consdata->nvars);
1826  assert(consdata->vars[pos] != NULL);
1827
1828  SCIP_CALL( SCIPdropVarEvent(scip, consdata->vars[pos],
1829  SCIP_EVENTTYPE_BOUNDTIGHTENED, eventhdlr, (SCIP_EVENTDATA*)consdata, -1) );
1830
1831  return SCIP_OKAY;
1832 }
1833
1834 /** drops bound change events for all variables in transformed linear constraint */
1835 static
1837  SCIP* scip, /**< SCIP data structure */
1838  SCIP_CONSDATA* consdata, /**< linear constraint data */
1839  SCIP_EVENTHDLR* eventhdlr /**< event handler to call for the event processing */
1840  )
1841 {
1842  int v;
1843
1844  assert(scip != NULL);
1845  assert(consdata != NULL);
1846
1847  /* drop event of every single variable */
1848  for( v = 0; v < consdata->nvars; ++v )
1849  {
1850  SCIP_CALL( consdataDropEvents(scip, consdata, eventhdlr, v) );
1851  }
1852
1853  return SCIP_OKAY;
1854 }
1855
1856 /** initialize variable lock data structure */
1857 static
1858 void initializeLocks(
1859  SCIP_CONSDATA* consdata, /**< constraint data */
1860  SCIP_Bool locked /**< should the variable be locked? */
1861  )
1862 {
1863  int nvars;
1864  int v;
1865
1866  nvars = consdata->nvars;
1867
1868  /* initialize locking arrays */
1869  for( v = 0; v < nvars; ++v )
1870  {
1871  consdata->downlocks[v] = locked;
1872  consdata->uplocks[v] = locked;
1873  }
1874 }
1875
1876 /** creates constraint data of cumulative constraint */
1877 static
1879  SCIP* scip, /**< SCIP data structure */
1880  SCIP_CONSDATA** consdata, /**< pointer to consdata */
1881  SCIP_VAR** vars, /**< array of integer variables */
1882  SCIP_CONS** linkingconss, /**< array of linking constraints for the integer variables, or NULL */
1883  int* durations, /**< array containing corresponding durations */
1884  int* demands, /**< array containing corresponding demands */
1885  int nvars, /**< number of variables */
1886  int capacity, /**< available cumulative capacity */
1887  int hmin, /**< left bound of time axis to be considered (including hmin) */
1888  int hmax, /**< right bound of time axis to be considered (not including hmax) */
1889  SCIP_Bool check /**< is the corresponding constraint a check constraint */
1890  )
1891 {
1892  int v;
1893
1894  assert(scip != NULL);
1895  assert(consdata != NULL);
1896  assert(vars != NULL || nvars > 0);
1897  assert(demands != NULL);
1898  assert(durations != NULL);
1899  assert(capacity >= 0);
1900  assert(hmin >= 0);
1901  assert(hmin < hmax);
1902
1903  /* create constraint data */
1904  SCIP_CALL( SCIPallocBlockMemory(scip, consdata) );
1905
1906  (*consdata)->hmin = hmin;
1907  (*consdata)->hmax = hmax;
1908
1909  (*consdata)->capacity = capacity;
1910  (*consdata)->demandrows = NULL;
1911  (*consdata)->demandrowssize = 0;
1912  (*consdata)->ndemandrows = 0;
1913  (*consdata)->scoverrows = NULL;
1914  (*consdata)->nscoverrows = 0;
1915  (*consdata)->scoverrowssize = 0;
1916  (*consdata)->bcoverrows = NULL;
1917  (*consdata)->nbcoverrows = 0;
1918  (*consdata)->bcoverrowssize = 0;
1919  (*consdata)->nvars = nvars;
1920  (*consdata)->varssize = nvars;
1921  (*consdata)->signature = 0;
1922  (*consdata)->validsignature = FALSE;
1923  (*consdata)->normalized = FALSE;
1924  (*consdata)->covercuts = FALSE;
1925  (*consdata)->propagated = FALSE;
1926  (*consdata)->varbounds = FALSE;
1927  (*consdata)->triedsolving = FALSE;
1928
1929  if( nvars > 0 )
1930  {
1931  assert(vars != NULL); /* for flexelint */
1932
1933  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->vars, vars, nvars) );
1934  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->demands, demands, nvars) );
1935  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->durations, durations, nvars) );
1936  (*consdata)->linkingconss = NULL;
1937
1938  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->downlocks, nvars) );
1939  SCIP_CALL( SCIPallocBlockMemoryArray(scip, &(*consdata)->uplocks, nvars) );
1940
1941  /* initialize variable lock data structure; the locks are only used if the contraint is a check constraint */
1942  initializeLocks(*consdata, check);
1943
1944  if( linkingconss != NULL )
1945  {
1946  SCIP_CALL( SCIPduplicateBlockMemoryArray(scip, &(*consdata)->linkingconss, linkingconss, nvars) );
1947  }
1948
1949  /* transform variables, if they are not yet transformed */
1950  if( SCIPisTransformed(scip) )
1951  {
1952  SCIPdebugMsg(scip, "get tranformed variables and constraints\n");
1953
1954  /* get transformed variables and do NOT captures these */
1955  SCIP_CALL( SCIPgetTransformedVars(scip, (*consdata)->nvars, (*consdata)->vars, (*consdata)->vars) );
1956
1957  /* multi-aggregated variables cannot be replaced by active variable; therefore we mark all variables for not
1958  * been multi-aggregated
1959  */
1960  for( v = 0; v < nvars; ++v )
1961  {
1962  SCIP_CALL( SCIPmarkDoNotMultaggrVar(scip, (*consdata)->vars[v]) );
1963  }
1964
1965  if( linkingconss != NULL )
1966  {
1967  /* get transformed constraints and captures these */
1968  SCIP_CALL( SCIPtransformConss(scip, (*consdata)->nvars, (*consdata)->linkingconss, (*consdata)->linkingconss) );
1969
1970  for( v = 0; v < nvars; ++v )
1972  }
1973  }
1974  }
1975  else
1976  {
1977  (*consdata)->vars = NULL;
1978  (*consdata)->downlocks = NULL;
1979  (*consdata)->uplocks = NULL;
1980  (*consdata)->demands = NULL;
1981  (*consdata)->durations = NULL;
1982  (*consdata)->linkingconss = NULL;
1983  }
1984
1985  /* initialize values for running propagation algorithms efficiently */
1986  (*consdata)->resstrength1 = -1.0;
1987  (*consdata)->resstrength2 = -1.0;
1988  (*consdata)->cumfactor1 = -1.0;
1989  (*consdata)->disjfactor1 = -1.0;
1990  (*consdata)->disjfactor2 = -1.0;
1991  (*consdata)->estimatedstrength = -1.0;
1992
1993  SCIPstatistic( (*consdata)->maxpeak = -1 );
1994
1995  return SCIP_OKAY;
1996 }
1997
1998 /** releases LP rows of constraint data and frees rows array */
1999 static
2001  SCIP* scip, /**< SCIP data structure */
2002  SCIP_CONSDATA** consdata /**< constraint data */
2003  )
2004 {
2005  int r;
2006
2007  assert(consdata != NULL);
2008  assert(*consdata != NULL);
2009
2010  for( r = 0; r < (*consdata)->ndemandrows; ++r )
2011  {
2012  assert((*consdata)->demandrows[r] != NULL);
2013  SCIP_CALL( SCIPreleaseRow(scip, &(*consdata)->demandrows[r]) );
2014  }
2015
2016  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->demandrows, (*consdata)->demandrowssize);
2017
2018  (*consdata)->ndemandrows = 0;
2019  (*consdata)->demandrowssize = 0;
2020
2021  /* free rows of cover cuts */
2022  for( r = 0; r < (*consdata)->nscoverrows; ++r )
2023  {
2024  assert((*consdata)->scoverrows[r] != NULL);
2025  SCIP_CALL( SCIPreleaseRow(scip, &(*consdata)->scoverrows[r]) );
2026  }
2027
2028  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->scoverrows, (*consdata)->scoverrowssize);
2029
2030  (*consdata)->nscoverrows = 0;
2031  (*consdata)->scoverrowssize = 0;
2032
2033  for( r = 0; r < (*consdata)->nbcoverrows; ++r )
2034  {
2035  assert((*consdata)->bcoverrows[r] != NULL);
2036  SCIP_CALL( SCIPreleaseRow(scip, &(*consdata)->bcoverrows[r]) );
2037  }
2038
2039  SCIPfreeBlockMemoryArrayNull(scip, &(*consdata)->bcoverrows, (*consdata)->bcoverrowssize);
2040
2041  (*consdata)->nbcoverrows = 0;
2042  (*consdata)->bcoverrowssize = 0;
2043
2044  (*consdata)->covercuts = FALSE;
2045
2046  return SCIP_OKAY;
2047 }
2048
2049 /** frees a cumulative constraint data */
2050 static
2052  SCIP* scip, /**< SCIP data structure */
2053  SCIP_CONSDATA** consdata /**< pointer to linear constraint data */
2054  )
2055 {
2057  int nvars;
2058
2059  assert(consdata != NULL);
2060  assert(*consdata != NULL);
2061
2062  nvars = (*consdata)->nvars;
2064
2065  if( varssize > 0 )
2066  {
2067  int v;
2068
2069  /* release and free the rows */
2070  SCIP_CALL( consdataFreeRows(scip, consdata) );
2071
2072  /* release the linking constraints if they were generated */
2073  if( (*consdata)->linkingconss != NULL )
2074  {
2075  for( v = nvars-1; v >= 0; --v )
2076  {
2077  assert((*consdata)->linkingconss[v] != NULL );
2078  SCIP_CALL( SCIPreleaseCons(scip, &(*consdata)->linkingconss[v]) );
2079  }
2080
2082  }
2083
2084  /* free arrays */
2085  SCIPfreeBlockMemoryArray(scip, &(*consdata)->downlocks, varssize);
2086  SCIPfreeBlockMemoryArray(scip, &(*consdata)->uplocks, varssize);
2087  SCIPfreeBlockMemoryArray(scip, &(*consdata)->durations, varssize);
2088  SCIPfreeBlockMemoryArray(scip, &(*consdata)->demands, varssize);
2089  SCIPfreeBlockMemoryArray(scip, &(*consdata)->vars, varssize);
2090  }
2091
2092  /* free memory */
2093  SCIPfreeBlockMemory(scip, consdata);
2094
2095  return SCIP_OKAY;
2096 }
2097
2098 /** prints cumulative constraint to file stream */
2099 static
2100 void consdataPrint(
2101  SCIP* scip, /**< SCIP data structure */
2102  SCIP_CONSDATA* consdata, /**< cumulative constraint data */
2103  FILE* file /**< output file (or NULL for standard output) */
2104  )
2105 {
2106  int v;
2107
2108  assert(consdata != NULL);
2109
2110  /* print coefficients */
2111  SCIPinfoMessage( scip, file, "cumulative(");
2112
2113  for( v = 0; v < consdata->nvars; ++v )
2114  {
2115  assert(consdata->vars[v] != NULL);
2116  if( v > 0 )
2117  SCIPinfoMessage(scip, file, ", ");
2118  SCIPinfoMessage(scip, file, "<%s>[%g,%g](%d)[%d]", SCIPvarGetName(consdata->vars[v]),
2119  SCIPvarGetLbGlobal(consdata->vars[v]), SCIPvarGetUbGlobal(consdata->vars[v]),
2120  consdata->durations[v], consdata->demands[v]);
2121  }
2122  SCIPinfoMessage(scip, file, ")[%d,%d) <= %d", consdata->hmin, consdata->hmax, consdata->capacity);
2123 }
2124
2125 /** deletes coefficient at given position from constraint data */
2126 static
2128  SCIP* scip, /**< SCIP data structure */
2129  SCIP_CONSDATA* consdata, /**< cumulative constraint data */
2130  SCIP_CONS* cons, /**< knapsack constraint */
2131  int pos /**< position of coefficient to delete */
2132  )
2133 {
2134  SCIP_CONSHDLR* conshdlr;
2135  SCIP_CONSHDLRDATA* conshdlrdata;
2136
2137  assert(scip != NULL);
2138  assert(consdata != NULL);
2139  assert(cons != NULL);
2140  assert(SCIPconsIsTransformed(cons));
2141  assert(!SCIPinProbing(scip));
2142
2143  SCIPdebugMsg(scip, "cumulative constraint <%s>: remove variable <%s>\n",
2144  SCIPconsGetName(cons), SCIPvarGetName(consdata->vars[pos]));
2145
2146  /* remove the rounding locks for the deleted variable */
2147  SCIP_CALL( SCIPunlockVarCons(scip, consdata->vars[pos], cons, consdata->downlocks[pos], consdata->uplocks[pos]) );
2148
2149  consdata->downlocks[pos] = FALSE;
2150  consdata->uplocks[pos] = FALSE;
2151
2152  if( consdata->linkingconss != NULL )
2153  {
2154  SCIP_CALL( SCIPreleaseCons(scip, &consdata->linkingconss[pos]) );
2155  }
2156
2157  /* get event handler */
2158  conshdlr = SCIPconsGetHdlr(cons);
2159  assert(conshdlr != NULL);
2160  conshdlrdata = SCIPconshdlrGetData(conshdlr);
2161  assert(conshdlrdata != NULL);
2162  assert(conshdlrdata->eventhdlr != NULL);
2163
2164  /* drop events */
2165  SCIP_CALL( consdataDropEvents(scip, consdata, conshdlrdata->eventhdlr, pos) );
2166
2167  SCIPdebugMsg(scip, "remove variable <%s>[%g,%g] from cumulative constraint <%s>\n",
2168  SCIPvarGetName(consdata->vars[pos]), SCIPvarGetLbGlobal(consdata->vars[pos]), SCIPvarGetUbGlobal(consdata->vars[pos]), SCIPconsGetName(cons));
2169
2170  /* in case the we did not remove the variable in the last slot of the arrays we move the current last to this
2171  * position
2172  */
2173  if( pos != consdata->nvars - 1 )
2174  {
2175  consdata->vars[pos] = consdata->vars[consdata->nvars-1];
2176  consdata->downlocks[pos] = consdata->downlocks[consdata->nvars-1];
2177  consdata->uplocks[pos] = consdata->uplocks[consdata->nvars-1];
2178  consdata->demands[pos] = consdata->demands[consdata->nvars-1];
2179  consdata->durations[pos] = consdata->durations[consdata->nvars-1];
2180
2181  if( consdata->linkingconss != NULL )
2182  {
2184  }
2185  }
2186
2187  consdata->nvars--;
2188  consdata->validsignature = FALSE;
2189  consdata->normalized = FALSE;
2190
2191  return SCIP_OKAY;
2192 }
2193
2194 /** collect linking constraints for each integer variable */
2195 static
2197  SCIP* scip, /**< SCIP data structure */
2198  SCIP_CONSDATA* consdata /**< pointer to consdata */
2199  )
2200 {
2201  int nvars;
2202  int v;
2203
2204  assert(scip != NULL);
2205  assert(consdata != NULL);
2206
2207  nvars = consdata->nvars;
2208  assert(nvars > 0);
2209  assert(consdata->linkingconss == NULL);
2210
2212
2213  for( v = 0; v < nvars; ++v )
2214  {
2215  SCIP_CONS* cons;
2216  SCIP_VAR* var;
2217
2218  var = consdata->vars[v];
2219  assert(var != NULL);
2220
2221  SCIPdebugMsg(scip, "linking constraint (%d of %d) for variable <%s>\n", v+1, nvars, SCIPvarGetName(var));
2222
2223  /* create linking constraint if it does not exist yet */
2224  if( !SCIPexistsConsLinking(scip, var) )
2225  {
2226  char name[SCIP_MAXSTRLEN];
2227
2228  (void)SCIPsnprintf(name, SCIP_MAXSTRLEN, "link(%s)", SCIPvarGetName(var));
2229
2230  /* creates and captures an linking constraint */
2231  SCIP_CALL( SCIPcreateConsLinking(scip, &cons, name, var, NULL, 0, 0,
2232  TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE /*TRUE*/, FALSE) );
2233  SCIP_CALL( SCIPaddCons(scip, cons) );
2234  consdata->linkingconss[v] = cons;
2235  }
2236  else
2237  {
2239  SCIP_CALL( SCIPcaptureCons(scip, consdata->linkingconss[v]) );
2240  }
2241
2243  assert(consdata->linkingconss[v] != NULL);
2246  }
2247
2248  return SCIP_OKAY;
2249 }
2250
2251 /**@} */
2252
2253
2254 /**@name Check methods
2255  *
2256  * @{
2257  */
2258
2259 /** check for the given starting time variables with their demands and durations if the cumulative conditions for the
2260  * given solution is satisfied
2261  */
2262 static
2264  SCIP* scip, /**< SCIP data structure */
2265  SCIP_SOL* sol, /**< primal solution, or NULL for current LP/pseudo solution */
2266  int nvars, /**< number of variables (jobs) */
2267  SCIP_VAR** vars, /**< array of integer variable which corresponds to starting times for a job */
2268  int* durations, /**< array containing corresponding durations */
2269  int* demands, /**< array containing corresponding demands */
2270  int capacity, /**< available cumulative capacity */
2271  int hmin, /**< left bound of time axis to be considered (including hmin) */
2272  int hmax, /**< right bound of time axis to be considered (not including hmax) */
2273  SCIP_Bool* violated, /**< pointer to store if the cumulative condition is violated */
2274  SCIP_CONS* cons, /**< constraint which is checked */
2275  SCIP_Bool printreason /**< should the reason for the violation be printed? */
2276  )
2277 {
2278  int* startsolvalues; /* stores when each job is starting */
2279  int* endsolvalues; /* stores when each job ends */
2280  int* startindices; /* we will sort the startsolvalues, thus we need to know which index of a job it corresponds to */
2281  int* endindices; /* we will sort the endsolvalues, thus we need to know which index of a job it corresponds to */
2282
2283  int freecapacity;
2284  int curtime; /* point in time which we are just checking */
2285  int endindex; /* index of endsolvalues with: endsolvalues[endindex] > curtime */
2286  int j;
2287
2288  SCIP_Real absviol;
2289  SCIP_Real relviol;
2290
2291  assert(scip != NULL);
2292  assert(violated != NULL);
2293
2294  (*violated) = FALSE;
2295
2296  if( nvars == 0 )
2297  return SCIP_OKAY;
2298
2299  assert(vars != NULL);
2300  assert(demands != NULL);
2301  assert(durations != NULL);
2302
2303  /* compute time points where we have to check whether capacity constraint is infeasible or not */
2304  SCIP_CALL( SCIPallocBufferArray(scip, &startsolvalues, nvars) );
2305  SCIP_CALL( SCIPallocBufferArray(scip, &endsolvalues, nvars) );
2306  SCIP_CALL( SCIPallocBufferArray(scip, &startindices, nvars) );
2307  SCIP_CALL( SCIPallocBufferArray(scip, &endindices, nvars) );
2308
2309  /* assign variables, start and endpoints to arrays */
2310  for ( j = 0; j < nvars; ++j )
2311  {
2312  int solvalue;
2313
2314  /* the constraint of the cumulative constraint handler should be called after the integrality check */
2315  assert(SCIPisFeasIntegral(scip, SCIPgetSolVal(scip, sol, vars[j])));
2316
2317  solvalue = SCIPconvertRealToInt(scip, SCIPgetSolVal(scip, sol, vars[j]));
2318
2319  /* we need to ensure that we check at least one time point during the effective horizon; therefore we project all
2320  * jobs which start before hmin to hmin
2321  */
2322  startsolvalues[j] = MAX(solvalue, hmin);
2323  startindices[j] = j;
2324
2325  endsolvalues[j] = MAX(solvalue + durations[j], hmin);
2326  endindices[j] = j;
2327  }
2328
2329  /* sort the arrays not-decreasing according to start solution values and end solution values (and sort the
2330  * corresponding indices in the same way)
2331  */
2332  SCIPsortIntInt(startsolvalues, startindices, nvars);
2333  SCIPsortIntInt(endsolvalues, endindices, nvars);
2334
2335  endindex = 0;
2336  freecapacity = capacity;
2337  absviol = 0.0;
2338  relviol = 0.0;
2339
2340  /* check each start point of a job whether the capacity is kept or not */
2341  for( j = 0; j < nvars; ++j )
2342  {
2343  /* only check intervals [hmin,hmax) */
2344  curtime = startsolvalues[j];
2345
2346  if( curtime >= hmax )
2347  break;
2348
2349  /* subtract all capacity needed up to this point */
2350  freecapacity -= demands[startindices[j]];
2351  while( j+1 < nvars && startsolvalues[j+1] == curtime )
2352  {
2353  j++;
2354  freecapacity -= demands[startindices[j]];
2355  }
2356
2357  /* free all capacity usages of jobs that are no longer running */
2358  while( endindex < nvars && curtime >= endsolvalues[endindex] )
2359  {
2360  freecapacity += demands[endindices[endindex]];
2361  ++endindex;
2362  }
2363  assert(freecapacity <= capacity);
2364
2365  /* update absolute and relative violation */
2366  if( absviol < (SCIP_Real) (-freecapacity) )
2367  {
2368  absviol = -freecapacity;
2369  relviol = SCIPrelDiff((SCIP_Real)(capacity - freecapacity), (SCIP_Real)capacity);
2370  }
2371
2372  /* check freecapacity to be smaller than zero */
2373  if( freecapacity < 0 && curtime >= hmin )
2374  {
2375  SCIPdebugMsg(scip, "freecapacity = %3d\n", freecapacity);
2376  (*violated) = TRUE;
2377
2378  if( printreason )
2379  {
2380  int i;
2381
2382  /* first state the violated constraints */
2383  SCIP_CALL( SCIPprintCons(scip, cons, NULL) );
2384
2385  /* second state the reason */
2386  SCIPinfoMessage(scip, NULL,
2387  ";\nviolation: at time point %d available capacity = %d, needed capacity = %d\n",
2388  curtime, capacity, capacity - freecapacity);
2389
2390  for( i = 0; i <= j; ++i )
2391  {
2392  if( startsolvalues[i] + durations[startindices[i]] > curtime )
2393  {
2394  SCIPinfoMessage(scip, NULL, "activity %s, start = %i, duration = %d, demand = %d \n",
2395  SCIPvarGetName(vars[startindices[i]]), startsolvalues[i], durations[startindices[i]],
2396  demands[startindices[i]]);
2397  }
2398  }
2399  }
2400  break;
2401  }
2402  } /*lint --e{850}*/
2403
2404  /* update constraint violation in solution */
2405  if( sol != NULL )
2406  SCIPupdateSolConsViolation(scip, sol, absviol, relviol);
2407
2408  /* free all buffer arrays */
2409  SCIPfreeBufferArray(scip, &endindices);
2410  SCIPfreeBufferArray(scip, &startindices);
2411  SCIPfreeBufferArray(scip, &endsolvalues);
2412  SCIPfreeBufferArray(scip, &startsolvalues);
2413
2414  return SCIP_OKAY;
2415 }
2416
2417 /** check if the given constrait is valid; checks each starting point of a job whether the remaining capacity is at
2418  * least zero or not. If not (*violated) is set to TRUE
2419  */
2420 static
2422  SCIP* scip, /**< SCIP data structure */
2423  SCIP_CONS* cons, /**< constraint to be checked */
2424  SCIP_SOL* sol, /**< primal solution, or NULL for current LP/pseudo solution */
2425  SCIP_Bool* violated, /**< pointer to store if the constraint is violated */
2426  SCIP_Bool printreason /**< should the reason for the violation be printed? */
2427  )
2428 {
2429  SCIP_CONSDATA* consdata;
2430
2431  assert(scip != NULL);
2432  assert(cons != NULL);
2433  assert(violated != NULL);
2434
2435  SCIPdebugMsg(scip, "check cumulative constraints <%s>\n", SCIPconsGetName(cons));
2436
2437  consdata = SCIPconsGetData(cons);
2438  assert(consdata != NULL);
2439
2440  /* check the cumulative condition */
2441  SCIP_CALL( checkCumulativeCondition(scip, sol, consdata->nvars, consdata->vars,
2442  consdata->durations, consdata->demands, consdata->capacity, consdata->hmin, consdata->hmax,
2443  violated, cons, printreason) );
2444
2445  return SCIP_OKAY;
2446 }
2447
2448 /**@} */
2449
2450 /**@name Conflict analysis
2451  *
2452  * @{
2453  */
2454
2455 /** resolves the propagation of the core time algorithm */
2456 static
2458  SCIP* scip, /**< SCIP data structure */
2459  int nvars, /**< number of start time variables (activities) */
2460  SCIP_VAR** vars, /**< array of start time variables */
2461  int* durations, /**< array of durations */
2462  int* demands, /**< array of demands */
2463  int capacity, /**< cumulative capacity */
2464  int hmin, /**< left bound of time axis to be considered (including hmin) */
2465  int hmax, /**< right bound of time axis to be considered (not including hmax) */
2466  SCIP_VAR* infervar, /**< inference variable */
2467  int inferdemand, /**< demand of the inference variable */
2468  int inferpeak, /**< time point which causes the propagation */
2469  int relaxedpeak, /**< relaxed time point which would be sufficient to be proved */
2470  SCIP_BDCHGIDX* bdchgidx, /**< the index of the bound change, representing the point of time where the change took place */
2471  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
2472  int* provedpeak, /**< pointer to store the actually proved peak, or NULL */
2473  SCIP_Bool* explanation /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
2474  )
2475 {
2476  SCIP_VAR* var;
2477  SCIP_Bool* reported;
2478  int duration;
2479  int maxlst;
2480  int minect;
2481  int ect;
2482  int lst;
2483  int j;
2484
2485  assert(SCIPgetStage(scip) == SCIP_STAGE_SOLVING || SCIPinProbing(scip));
2486
2487  SCIPdebugMsg(scip, "variable <%s>: (demand %d) resolve propagation of core time algorithm (peak %d)\n",
2488  SCIPvarGetName(infervar), inferdemand, inferpeak);
2489  assert(nvars > 0);
2490
2491  /* adjusted capacity */
2492  capacity -= inferdemand;
2493  maxlst = INT_MIN;
2494  minect = INT_MAX;
2495
2496  SCIP_CALL( SCIPallocBufferArray(scip, &reported, nvars) );
2497  BMSclearMemoryArray(reported, nvars);
2498
2499  /* first we loop over all variables and adjust the capacity with those jobs which provide a global core at the
2500  * inference peak and those where the current conflict bounds provide a core at the inference peak
2501  */
2502  for( j = 0; j < nvars && capacity >= 0; ++j )
2503  {
2504  var = vars[j];
2505  assert(var != NULL);
2506
2507  /* skip inference variable */
2508  if( var == infervar )
2509  continue;
2510
2511  duration = durations[j];
2512  assert(duration > 0);
2513
2514  /* compute cores of jobs; if core overlaps interval of inference variable add this job to the array */
2515  assert(!SCIPvarIsActive(var) || SCIPisFeasEQ(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, TRUE), SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE)));
2516  assert(SCIPisFeasIntegral(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, TRUE)));
2517  assert(!SCIPvarIsActive(var) || SCIPisFeasEQ(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, TRUE), SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)));
2518  assert(SCIPisFeasIntegral(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, TRUE)));
2519
2520  SCIPdebugMsg(scip, "variable <%s>: glb=[%g,%g] conflict=[%g,%g] (duration %d, demand %d)\n",
2522  SCIPgetConflictVarLb(scip, var), SCIPgetConflictVarUb(scip, var), duration, demands[j]);
2523
2524  ect = SCIPconvertRealToInt(scip, SCIPvarGetLbGlobal(var)) + duration;
2525  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbGlobal(var));
2526
2527  /* check if the inference peak is part of the global bound core; if so we decreasing the capacity by the demand of
2528  * that job without adding it the explanation
2529  */
2530  if( inferpeak < ect && lst <= inferpeak )
2531  {
2532  capacity -= demands[j];
2533  reported[j] = TRUE;
2534
2535  maxlst = MAX(maxlst, lst);
2536  minect = MIN(minect, ect);
2537  assert(maxlst < minect);
2538
2539  if( explanation != NULL )
2540  explanation[j] = TRUE;
2541
2542  continue;
2543  }
2544
2545  /* collect the conflict bound core (the conflict bounds are those bounds which are already part of the conflict)
2546  * hence these bound are already reported by other resolve propation steps. In case a bound (lower or upper) is
2547  * not part of the conflict yet we get the global bounds back.
2548  */
2549  ect = SCIPconvertRealToInt(scip, SCIPgetConflictVarLb(scip, var)) + duration;
2550  lst = SCIPconvertRealToInt(scip, SCIPgetConflictVarUb(scip, var));
2551
2552  /* check if the inference peak is part of the conflict bound core; if so we decreasing the capacity by the demand
2553  * of that job without and collect the job as part of the explanation
2554  *
2555  * @note we do not need to reported that job to SCIP since the required bounds are already reported
2556  */
2557  if( inferpeak < ect && lst <= inferpeak )
2558  {
2559  capacity -= demands[j];
2560  reported[j] = TRUE;
2561
2562  maxlst = MAX(maxlst, lst);
2563  minect = MIN(minect, ect);
2564  assert(maxlst < minect);
2565
2566  if( explanation != NULL )
2567  explanation[j] = TRUE;
2568  }
2569  }
2570
2571  if( capacity >= 0 )
2572  {
2573  int* cands;
2574  int* canddemands;
2575  int ncands;
2576  int c;
2577
2578  SCIP_CALL( SCIPallocBufferArray(scip, &cands, nvars) );
2579  SCIP_CALL( SCIPallocBufferArray(scip, &canddemands, nvars) );
2580  ncands = 0;
2581
2582  /* collect all cores of the variables which lay in the considered time window except the inference variable */
2583  for( j = 0; j < nvars; ++j )
2584  {
2585  var = vars[j];
2586  assert(var != NULL);
2587
2588  /* skip inference variable */
2589  if( var == infervar || reported[j] )
2590  continue;
2591
2592  duration = durations[j];
2593  assert(duration > 0);
2594
2595  /* compute cores of jobs; if core overlaps interval of inference variable add this job to the array */
2596  assert(!SCIPvarIsActive(var) || SCIPisFeasEQ(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, TRUE), SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE)));
2597  assert(SCIPisFeasIntegral(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, TRUE)));
2598  assert(!SCIPvarIsActive(var) || SCIPisFeasEQ(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, TRUE), SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)));
2599  assert(SCIPisFeasIntegral(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, TRUE)));
2600
2601  /* collect local core information */
2602  ect = SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)) + duration;
2603  lst = SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE));
2604
2605  SCIPdebugMsg(scip, "variable <%s>: loc=[%g,%g] glb=[%g,%g] (duration %d, demand %d)\n",
2606  SCIPvarGetName(var), SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE), SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE),
2607  SCIPvarGetLbGlobal(var), SCIPvarGetUbGlobal(var), duration, demands[j]);
2608
2609  /* check if the inference peak is part of the core */
2610  if( inferpeak < ect && lst <= inferpeak )
2611  {
2612  cands[ncands] = j;
2613  canddemands[ncands] = demands[j];
2614  ncands++;
2615
2616  capacity -= demands[j];
2617  }
2618  }
2619
2620  /* sort candidates indices w.r.t. their demands */
2621  SCIPsortDownIntInt(canddemands, cands, ncands);
2622
2623  assert(capacity < 0);
2624  assert(ncands > 0);
2625
2626  /* greedily remove candidates form the list such that the needed capacity is still exceeded */
2627  while( capacity + canddemands[ncands-1] < 0 )
2628  {
2629  ncands--;
2630  capacity += canddemands[ncands];
2631  assert(ncands > 0);
2632  }
2633
2634  /* compute the size (number of time steps) of the job cores */
2635  for( c = 0; c < ncands; ++c )
2636  {
2637  var = vars[cands[c]];
2638  assert(var != NULL);
2639
2640  duration = durations[cands[c]];
2641
2642  ect = SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)) + duration;
2643  lst = SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE));
2644
2645  maxlst = MAX(maxlst, lst);
2646  minect = MIN(minect, ect);
2647  assert(maxlst < minect);
2648  }
2649
2650  SCIPdebugMsg(scip, "infer peak %d, relaxed peak %d, lst %d, ect %d\n", inferpeak, relaxedpeak, maxlst, minect);
2651  assert(inferpeak >= maxlst);
2652  assert(inferpeak < minect);
2653
2654  /* check if the collect variable are sufficient to prove the relaxed bound (relaxedpeak) */
2655  if( relaxedpeak < inferpeak )
2656  {
2657  inferpeak = MAX(maxlst, relaxedpeak);
2658  }
2659  else if( relaxedpeak > inferpeak )
2660  {
2661  inferpeak = MIN(minect-1, relaxedpeak);
2662  }
2663  assert(inferpeak >= hmin);
2664  assert(inferpeak < hmax);
2665  assert(inferpeak >= maxlst);
2666  assert(inferpeak < minect);
2667
2668  /* post all necessary bound changes */
2669  for( c = 0; c < ncands; ++c )
2670  {
2671  var = vars[cands[c]];
2672  assert(var != NULL);
2673
2674  if( usebdwidening )
2675  {
2676  duration = durations[cands[c]];
2677
2678  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, var, bdchgidx, (SCIP_Real)(inferpeak - duration + 1)) );
2679  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, var, bdchgidx, (SCIP_Real)inferpeak) );
2680  }
2681  else
2682  {
2683  SCIP_CALL( SCIPaddConflictLb(scip, var, bdchgidx) );
2684  SCIP_CALL( SCIPaddConflictUb(scip, var, bdchgidx) );
2685  }
2686
2687  if( explanation != NULL )
2688  explanation[cands[c]] = TRUE;
2689  }
2690
2691  SCIPfreeBufferArray(scip, &canddemands);
2692  SCIPfreeBufferArray(scip, &cands);
2693  }
2694
2695  SCIPfreeBufferArray(scip, &reported);
2696
2697  if( provedpeak != NULL )
2698  *provedpeak = inferpeak;
2699
2700  return SCIP_OKAY;
2701 }
2702
2703 #if 0
2704 /** repropagation of edge finding algorithm simplified version from Petr Vilim only a small subset is reported such that
2705  * energy in total and for bound change is enough
2706  */
2707 static
2708 SCIP_RETCODE resolvePropagationEdgeFinding(
2709  SCIP* scip, /**< SCIP data structure */
2710  int nvars, /**< number of start time variables (activities) */
2711  SCIP_VAR** vars, /**< array of start time variables */
2712  int* durations, /**< array of durations */
2713  int hmin, /**< left bound of time axis to be considered (including hmin) */
2714  int hmax, /**< right bound of time axis to be considered (not including hmax) */
2715  SCIP_VAR* infervar, /**< variable whose bound change is to be explained */
2716  INFERINFO inferinfo, /**< inference info containing position of correct bdchgids */
2717  SCIP_BOUNDTYPE boundtype, /**< the type of the changed bound (lower or upper bound) */
2718  SCIP_BDCHGIDX* bdchgidx, /**< the index of the bound change, representing the point of time where the change took place */
2719  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
2720  SCIP_Bool* explanation /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
2721  )
2722 {
2723  int est;
2724  int lct;
2725  int j;
2726
2727  /* ???????????????????? do bound widening */
2728
2729  assert(scip != NULL);
2730  assert(nvars > 0);
2731  assert(inferInfoGetProprule(inferinfo) == PROPRULE_2_EDGEFINDING);
2732
2733  SCIPdebugMsg(scip, "repropagate edge-finding with short reasons for variable <%s>\n", SCIPvarGetName(infervar));
2734
2735  if( boundtype == SCIP_BOUNDTYPE_LOWER )
2736  {
2737  SCIP_CALL( SCIPaddConflictLb(scip, infervar, bdchgidx) );
2738  }
2739  else
2740  {
2741  SCIP_CALL( SCIPaddConflictUb(scip, infervar, bdchgidx) );
2742  }
2743
2744  est = inferInfoGetData1(inferinfo);
2745  lct = inferInfoGetData2(inferinfo);
2746  assert(est < lct);
2747
2748  /* collect the energies of all variables in [est_omega, lct_omega] */
2749  for( j = 0; j < nvars; ++j )
2750  {
2751  SCIP_VAR* var;
2752  SCIP_Bool left;
2753  SCIP_Bool right;
2754  int duration;
2755  int lb;
2756  int ub;
2757
2758  var = vars[j];
2759  assert(var != NULL);
2760
2761  if( var == infervar )
2762  {
2763  if( explanation != NULL )
2764  explanation[j] = TRUE;
2765
2766  continue;
2767  }
2768
2769  lb = SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE));
2770  ub = SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE));
2771
2772  duration = durations[j];
2773  assert(duration > 0);
2774
2775  /* in case the earliest start time is equal to hmin we have to also consider the jobs which run in that region
2776  * since we use adjusted jobs during the propagation
2777  */
2778  left = (est == hmin && lb + duration > hmin) || lb >= est;
2779
2780  /* in case the latest completion time is equal to hmax we have to also consider the jobs which run in that region
2781  * since we use adjusted jobs during the propagation
2782  */
2783  right = (lct == hmax && ub < hmax) || ub + duration <= lct;
2784
2785  /* store all jobs running in [est_omega; lct_omega] */
2786  if( left && right )
2787  {
2788  /* check if bound widening should be used */
2789  if( usebdwidening )
2790  {
2791  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, var, bdchgidx, (SCIP_Real)(lct - duration)) );
2792  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, var, bdchgidx, (SCIP_Real)(est)) );
2793  }
2794  else
2795  {
2796  SCIP_CALL( SCIPaddConflictLb(scip, var, bdchgidx) );
2797  SCIP_CALL( SCIPaddConflictUb(scip, var, bdchgidx) );
2798  }
2799
2800  if( explanation != NULL )
2801  explanation[j] = TRUE;
2802  }
2803  }
2804
2805  return SCIP_OKAY;
2806 }
2807 #endif
2808
2809 /** compute the minimum overlaps w.r.t. the duration of the job and the time window [begin,end) */
2810 static
2811 int computeOverlap(
2812  int begin, /**< begin of the times interval */
2813  int end, /**< end of time interval */
2814  int est, /**< earliest start time */
2815  int lst, /**< latest start time */
2816  int duration /**< duration of the job */
2817  )
2818 {
2819  int left;
2820  int right;
2821  int ect;
2822  int lct;
2823
2824  ect = est + duration;
2825  lct = lst + duration;
2826
2827  /* check if job runs completely within [begin,end) */
2828  if( lct <= end && est >= begin )
2829  return duration;
2830
2831  assert(lst <= end && ect >= begin);
2832
2833  left = ect - begin;
2834  assert(left > 0);
2835
2836  right = end - lst;
2837  assert(right > 0);
2838
2839  return MIN3(left, right, end - begin);
2840 }
2841
2842 /** an overload was detected due to the time-time edge-finding propagate; initialized conflict analysis, add an initial
2843  * reason
2844  *
2845  * @note the conflict analysis is not performend, only the initialized SCIP_Bool pointer is set to TRUE
2846  */
2847 static
2849  SCIP* scip, /**< SCIP data structure */
2850  int nvars, /**< number of start time variables (activities) */
2851  SCIP_VAR** vars, /**< array of start time variables */
2852  int* durations, /**< array of durations */
2853  int* demands, /**< array of demands */
2854  int capacity, /**< capacity of the cumulative condition */
2855  int begin, /**< begin of the time window */
2856  int end, /**< end of the time window */
2857  SCIP_VAR* infervar, /**< variable which was propagate, or NULL */
2858  SCIP_BOUNDTYPE boundtype, /**< the type of the changed bound (lower or upper bound) */
2859  SCIP_BDCHGIDX* bdchgidx, /**< the index of the bound change, representing the point of time where the change took place */
2860  SCIP_Real relaxedbd, /**< the relaxed bound which is sufficient to be explained */
2861  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
2862  SCIP_Bool* explanation /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
2863  )
2864 {
2865  int* locenergies;
2866  int* overlaps;
2867  int* idxs;
2868
2869  SCIP_Longint requiredenergy;
2870  int v;
2871
2872  SCIP_CALL( SCIPallocBufferArray(scip, &locenergies, nvars) );
2873  SCIP_CALL( SCIPallocBufferArray(scip, &overlaps, nvars) );
2874  SCIP_CALL( SCIPallocBufferArray(scip, &idxs, nvars) );
2875
2876  /* energy which needs be explained */
2877  requiredenergy = ((SCIP_Longint) end - begin) * capacity;
2878
2879  SCIPdebugMsg(scip, "analysis energy load in [%d,%d) (capacity %d, energy %d)\n", begin, end, capacity, requiredenergy);
2880
2881  /* collect global contribution and adjusted the required energy by the amount of energy the inference variable
2882  * takes
2883  */
2884  for( v = 0; v < nvars; ++v )
2885  {
2886  SCIP_VAR* var;
2887  int glbenergy;
2888  int duration;
2889  int demand;
2890  int est;
2891  int lst;
2892
2893  var = vars[v];
2894  assert(var != NULL);
2895
2896  locenergies[v] = 0;
2897  overlaps[v] = 0;
2898  idxs[v] = v;
2899
2900  demand = demands[v];
2901  assert(demand > 0);
2902
2903  duration = durations[v];
2904  assert(duration > 0);
2905
2906  /* check if the variable equals the inference variable (the one which was propagated) */
2907  if( infervar == var )
2908  {
2909  int overlap;
2910  int right;
2911  int left;
2912
2913  assert(relaxedbd != SCIP_UNKNOWN); /*lint !e777*/
2914
2915  SCIPdebugMsg(scip, "inference variable <%s>[%g,%g] %s %g (duration %d, demand %d)\n",
2916  SCIPvarGetName(var), SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE), SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE),
2917  boundtype == SCIP_BOUNDTYPE_LOWER ? ">=" : "<=", relaxedbd, duration, demand);
2918
2919  /* compute the amount of energy which needs to be available for enforcing the propagation and report the bound
2920  * which is necessary from the inference variable
2921  */
2922  if( boundtype == SCIP_BOUNDTYPE_UPPER )
2923  {
2924  int lct;
2925
2926  /* get the latest start time of the infer start time variable before the propagation took place */
2927  lst = SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE));
2928
2929  /* the latest start time of the inference start time variable before the propagation needs to be smaller as
2930  * the end of the time interval; meaning the job needs be overlap with the time interval in case the job is
2931  * scheduled w.r.t. its latest start time
2932  */
2933  assert(lst < end);
2934
2935  /* compute the overlap of the job in case it would be scheduled w.r.t. its latest start time and the time
2936  * interval (before the propagation)
2937  */
2938  right = MIN3(end - lst, end - begin, duration);
2939
2940  /* the job needs to overlap with the interval; otherwise the propagation w.r.t. this time window is not valid */
2941  assert(right > 0);
2942
2943  lct = SCIPconvertRealToInt(scip, relaxedbd) + duration;
2944  assert(begin <= lct);
2945  assert(bdchgidx == NULL || SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, TRUE)) < begin);
2946
2947  /* compute the overlap of the job after the propagation but considering the relaxed bound */
2948  left = MIN(lct - begin + 1, end - begin);
2949  assert(left > 0);
2950
2951  /* compute the minimum overlap; */
2952  overlap = MIN(left, right);
2953  assert(overlap > 0);
2954  assert(overlap <= end - begin);
2955  assert(overlap <= duration);
2956
2957  if( usebdwidening )
2958  {
2959  assert(SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE)) <= (end - overlap));
2960  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, var, bdchgidx, (SCIP_Real)(end - overlap)) );
2961  }
2962  else
2963  {
2964  SCIP_CALL( SCIPaddConflictUb(scip, var, bdchgidx) );
2965  }
2966  }
2967  else
2968  {
2969  int ect;
2970
2971  assert(boundtype == SCIP_BOUNDTYPE_LOWER);
2972
2973  /* get the earliest completion time of the infer start time variable before the propagation took place */
2974  ect = SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)) + duration;
2975
2976  /* the earliest start time of the inference start time variable before the propagation needs to be larger as
2977  * than the beginning of the time interval; meaning the job needs be overlap with the time interval in case
2978  * the job is scheduled w.r.t. its earliest start time
2979  */
2980  assert(ect > begin);
2981
2982  /* compute the overlap of the job in case it would be scheduled w.r.t. its earliest start time and the time
2983  * interval (before the propagation)
2984  */
2985  left = MIN3(ect - begin, end - begin, duration);
2986
2987  /* the job needs to overlap with the interval; otherwise the propagation w.r.t. this time window is not valid */
2988  assert(left > 0);
2989
2990  est = SCIPconvertRealToInt(scip, relaxedbd);
2991  assert(end >= est);
2992  assert(bdchgidx == NULL || end - SCIPgetVarLbAtIndex(scip, var, bdchgidx, TRUE) < duration);
2993
2994  /* compute the overlap of the job after the propagation but considering the relaxed bound */
2995  right = MIN(end - est + 1, end - begin);
2996  assert(right > 0);
2997
2998  /* compute the minimum overlap */
2999  overlap = MIN(left, right);
3000  assert(overlap > 0);
3001  assert(overlap <= end - begin);
3002  assert(overlap <= duration);
3003
3004  if( usebdwidening )
3005  {
3006  assert(SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE)) >= (begin + overlap - duration));
3007  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, var, bdchgidx, (SCIP_Real)(begin + overlap - duration)) );
3008  }
3009  else
3010  {
3011  SCIP_CALL( SCIPaddConflictLb(scip, var, bdchgidx) );
3012  }
3013  }
3014
3015  /* subtract the amount of energy which is available due to the overlap of the inference start time */
3016  requiredenergy -= (SCIP_Longint) overlap * demand;
3017
3018  if( explanation != NULL )
3019  explanation[v] = TRUE;
3020
3021  continue;
3022  }
3023
3024  /* global time points */
3025  est = SCIPconvertRealToInt(scip, SCIPvarGetLbGlobal(var));
3026  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbGlobal(var));
3027
3028  glbenergy = 0;
3029
3030  /* check if the has any overlap w.r.t. global bound; meaning some parts of the job will run for sure within the
3031  * time window
3032  */
3033  if( est + duration > begin && lst < end )
3034  {
3035  /* evaluated global contribution */
3036  glbenergy = computeOverlap(begin, end, est, lst, duration) * demand;
3037
3038  /* remove the globally available energy form the required energy */
3039  requiredenergy -= glbenergy;
3040
3041  if( explanation != NULL )
3042  explanation[v] = TRUE;
3043  }
3044
3045  /* local time points */
3046  est = SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, var, bdchgidx, FALSE));
3047  lst = SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, var, bdchgidx, FALSE));
3048
3049  /* check if the job has any overlap w.r.t. local bound; meaning some parts of the job will run for sure within the
3050  * time window
3051  */
3052  if( est + duration > begin && lst < end )
3053  {
3054  overlaps[v] = computeOverlap(begin, end, est, lst, duration);
3055
3056  /* evaluated additionally local energy contribution */
3057  locenergies[v] = overlaps[v] * demand - glbenergy;
3058  assert(locenergies[v] >= 0);
3059  }
3060  }
3061
3062  /* sort the variable contributions w.r.t. additional local energy contributions */
3063  SCIPsortDownIntIntInt(locenergies, overlaps, idxs, nvars);
3064
3065  /* add local energy contributions until an overload is implied */
3066  for( v = 0; v < nvars && requiredenergy >= 0; ++v )
3067  {
3068  SCIP_VAR* var;
3069  int duration;
3070  int overlap;
3071  int relaxlb;
3072  int relaxub;
3073  int idx;
3074
3075  idx = idxs[v];
3076  assert(idx >= 0 && idx < nvars);
3077
3078  var = vars[idx];
3079  assert(var != NULL);
3080  assert(var != infervar);
3081
3082  duration = durations[idx];
3083  assert(duration > 0);
3084
3085  overlap = overlaps[v];
3086  assert(overlap > 0);
3087
3088  requiredenergy -= locenergies[v];
3089
3090  if( requiredenergy < -1 )
3091  {
3092  int demand;
3093
3094  demand = demands[idx];
3095  assert(demand > 0);
3096
3097  overlap += (int)((requiredenergy + 1) / demand);
3098
3099 #ifndef NDEBUG
3100  requiredenergy += locenergies[v];
3101  requiredenergy -= (SCIP_Longint) overlap * demand;
3102  assert(requiredenergy < 0);
3103 #endif
3104  }
3105  assert(overlap > 0);
3106
3107  relaxlb = begin - duration + overlap;
3108  relaxub = end - overlap;
3109
3110  SCIPdebugMsg(scip, "variable <%s> glb=[%g,%g] loc=[%g,%g], conf=[%g,%g], added=[%d,%d] (demand %d, duration %d)\n",
3111  SCIPvarGetName(var),
3114  SCIPgetConflictVarLb(scip, var), SCIPgetConflictVarUb(scip, var),
3115  relaxlb, relaxub, demands[idx], duration);
3116
3117  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, var, bdchgidx, (SCIP_Real)relaxlb) );
3118  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, var, bdchgidx, (SCIP_Real)relaxub) );
3119
3120  if( explanation != NULL )
3121  explanation[idx] = TRUE;
3122  }
3123
3124  assert(requiredenergy < 0);
3125
3126  SCIPfreeBufferArray(scip, &idxs);
3127  SCIPfreeBufferArray(scip, &overlaps);
3128  SCIPfreeBufferArray(scip, &locenergies);
3129
3130  return SCIP_OKAY;
3131 }
3132
3133 /** resolve propagation w.r.t. the cumulative condition */
3134 static
3136  SCIP* scip, /**< SCIP data structure */
3137  int nvars, /**< number of start time variables (activities) */
3138  SCIP_VAR** vars, /**< array of start time variables */
3139  int* durations, /**< array of durations */
3140  int* demands, /**< array of demands */
3141  int capacity, /**< cumulative capacity */
3142  int hmin, /**< left bound of time axis to be considered (including hmin) */
3143  int hmax, /**< right bound of time axis to be considered (not including hmax) */
3144  SCIP_VAR* infervar, /**< the conflict variable whose bound change has to be resolved */
3145  INFERINFO inferinfo, /**< the user information */
3146  SCIP_BOUNDTYPE boundtype, /**< the type of the changed bound (lower or upper bound) */
3147  SCIP_BDCHGIDX* bdchgidx, /**< the index of the bound change, representing the point of time where the change took place */
3148  SCIP_Real relaxedbd, /**< the relaxed bound which is sufficient to be explained */
3149  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
3150  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
3151  SCIP_RESULT* result /**< pointer to store the result of the propagation conflict resolving call */
3152  )
3153 {
3154  switch( inferInfoGetProprule(inferinfo) )
3155  {
3156  case PROPRULE_1_CORETIMES:
3157  {
3158  int inferdemand;
3159  int inferduration;
3160  int inferpos;
3161  int inferpeak;
3162  int relaxedpeak;
3163  int provedpeak;
3164
3165  /* get the position of the inferred variable in the vars array */
3166  inferpos = inferInfoGetData1(inferinfo);
3167  if( inferpos >= nvars || vars[inferpos] != infervar )
3168  {
3169  /* find inference variable in constraint */
3170  for( inferpos = 0; inferpos < nvars && vars[inferpos] != infervar; ++inferpos )
3171  {}
3172  }
3173  assert(inferpos < nvars);
3174  assert(vars[inferpos] == infervar);
3175
3176  inferdemand = demands[inferpos];
3177  inferduration = durations[inferpos];
3178
3179  if( boundtype == SCIP_BOUNDTYPE_UPPER )
3180  {
3181  /* we propagated the latest start time (upper bound) step wise with a step length of at most the duration of
3182  * the inference variable
3183  */
3184  assert(SCIPgetVarUbAtIndex(scip, infervar, bdchgidx, FALSE) - SCIPgetVarUbAtIndex(scip, infervar, bdchgidx, TRUE) < inferduration + 0.5);
3185
3186  SCIPdebugMsg(scip, "variable <%s>: upper bound changed from %g to %g (relaxed %g)\n",
3187  SCIPvarGetName(infervar), SCIPgetVarUbAtIndex(scip, infervar, bdchgidx, FALSE),
3188  SCIPgetVarUbAtIndex(scip, infervar, bdchgidx, TRUE), relaxedbd);
3189
3190  /* get the inference peak that the time point which lead to the that propagtion */
3191  inferpeak = inferInfoGetData2(inferinfo);
3192  /* the bound passed back to be resolved might be tighter as the bound propagted by the core time propagator;
3193  * this can happen if the variable is not activ and aggregated to an activ variable with a scale != 1.0
3194  */
3195  assert(SCIPconvertRealToInt(scip, SCIPgetVarUbAtIndex(scip, infervar, bdchgidx, TRUE)) + inferduration <= inferpeak);
3196  relaxedpeak = SCIPconvertRealToInt(scip, relaxedbd) + inferduration;
3197
3198  /* make sure that the relaxed peak is part of the effective horizon */
3199  relaxedpeak = MIN(relaxedpeak, hmax-1);
3200
3201  /* make sure that relaxed peak is not larger than the infer peak
3202  *
3203  * This can happen in case the variable is not an active variable!
3204  */
3205  relaxedpeak = MAX(relaxedpeak, inferpeak);
3206  assert(relaxedpeak >= inferpeak);
3207  assert(relaxedpeak >= hmin);
3208  }
3209  else
3210  {
3211  assert(boundtype == SCIP_BOUNDTYPE_LOWER);
3212
3213  SCIPdebugMsg(scip, "variable <%s>: lower bound changed from %g to %g (relaxed %g)\n",
3214  SCIPvarGetName(infervar), SCIPgetVarLbAtIndex(scip, infervar, bdchgidx, FALSE),
3215  SCIPgetVarLbAtIndex(scip, infervar, bdchgidx, TRUE), relaxedbd);
3216
3217  /* get the time interval where the job could not be scheduled */
3218  inferpeak = inferInfoGetData2(inferinfo);
3219  /* the bound passed back to be resolved might be tighter as the bound propagted by the core time propagator;
3220  * this can happen if the variable is not activ and aggregated to an activ variable with a scale != 1.0
3221  */
3222  assert(SCIPconvertRealToInt(scip, SCIPgetVarLbAtIndex(scip, infervar, bdchgidx, TRUE)) - 1 >= inferpeak);
3223  relaxedpeak = SCIPconvertRealToInt(scip, relaxedbd) - 1;
3224
3225  /* make sure that the relaxed peak is part of the effective horizon */
3226  relaxedpeak = MAX(relaxedpeak, hmin);
3227
3228  /* make sure that relaxed peak is not larger than the infer peak
3229  *
3230  * This can happen in case the variable is not an active variable!
3231  */
3232  relaxedpeak = MIN(relaxedpeak, inferpeak);
3233  assert(relaxedpeak < hmax);
3234  }
3235
3236  /* resolves the propagation of the core time algorithm */
3237  SCIP_CALL( resolvePropagationCoretimes(scip, nvars, vars, durations, demands, capacity, hmin, hmax,
3238  infervar, inferdemand, inferpeak, relaxedpeak, bdchgidx, usebdwidening, &provedpeak, explanation) );
3239
3240  if( boundtype == SCIP_BOUNDTYPE_UPPER )
3241  {
3242  if( usebdwidening )
3243  {
3244  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, infervar, NULL, (SCIP_Real)provedpeak) );
3245  }
3246  else
3247  {
3248  /* old upper bound of variable itself is part of the explanation */
3249  SCIP_CALL( SCIPaddConflictUb(scip, infervar, bdchgidx) );
3250  }
3251  }
3252  else
3253  {
3254  assert(boundtype == SCIP_BOUNDTYPE_LOWER);
3255
3256  if( usebdwidening )
3257  {
3258  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, infervar, bdchgidx, (SCIP_Real)(provedpeak - inferduration + 1)) );
3259  }
3260  else
3261  {
3262  /* old lower bound of variable itself is part of the explanation */
3263  SCIP_CALL( SCIPaddConflictLb(scip, infervar, bdchgidx) );
3264  }
3265  }
3266
3267  if( explanation != NULL )
3268  explanation[inferpos] = TRUE;
3269
3270  break;
3271  }
3273  case PROPRULE_3_TTEF:
3274  {
3275  int begin;
3276  int end;
3277
3278  begin = inferInfoGetData1(inferinfo);
3279  end = inferInfoGetData2(inferinfo);
3280  assert(begin < end);
3281
3282  begin = MAX(begin, hmin);
3283  end = MIN(end, hmax);
3284
3285  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
3286  begin, end, infervar, boundtype, bdchgidx, relaxedbd, usebdwidening, explanation) );
3287
3288  break;
3289  }
3290
3291  default:
3292  SCIPerrorMessage("invalid inference information %d\n", inferInfoGetProprule(inferinfo));
3293  SCIPABORT();
3294  return SCIP_INVALIDDATA; /*lint !e527*/
3295  }
3296
3297  (*result) = SCIP_SUCCESS;
3298
3299  return SCIP_OKAY;
3300 }
3301
3302 /**@} */
3303
3304
3305 /**@name Enforcement methods
3306  *
3307  * @{
3308  */
3309
3310 /** apply all fixings which are given by the alternative bounds */
3311 static
3313  SCIP* scip, /**< SCIP data structure */
3314  SCIP_VAR** vars, /**< array of active variables */
3315  int nvars, /**< number of active variables */
3316  int* alternativelbs, /**< alternative lower bounds */
3317  int* alternativeubs, /**< alternative lower bounds */
3318  int* downlocks, /**< number of constraints with down lock participating by the computation */
3319  int* uplocks, /**< number of constraints with up lock participating by the computation */
3320  SCIP_Bool* branched /**< pointer to store if a branching was applied */
3321  )
3322 {
3323  int v;
3324
3325  for( v = 0; v < nvars; ++v )
3326  {
3327  SCIP_VAR* var;
3328  SCIP_Real objval;
3329
3330  var = vars[v];
3331  assert(var != NULL);
3332
3333  objval = SCIPvarGetObj(var);
3334
3335  if( SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL) == downlocks[v] && !SCIPisNegative(scip, objval) )
3336  {
3337  int ub;
3338
3339  ub = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
3340
3341  if( alternativelbs[v] <= ub )
3342  {
3343  SCIP_CALL( SCIPbranchVarHole(scip, var, SCIPvarGetLbLocal(var), (SCIP_Real)alternativelbs[v], NULL, NULL) );
3344  (*branched) = TRUE;
3345
3346  SCIPdebugMsg(scip, "variable <%s> branched domain hole (%g,%d)\n", SCIPvarGetName(var),
3347  SCIPvarGetLbLocal(var), alternativelbs[v]);
3348
3349  return SCIP_OKAY;
3350  }
3351  }
3352
3353  if( SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) == uplocks[v] && !SCIPisPositive(scip, objval) )
3354  {
3355  int lb;
3356
3357  lb = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
3358
3359  if( alternativeubs[v] >= lb )
3360  {
3361  SCIP_CALL( SCIPbranchVarHole(scip, var, (SCIP_Real)alternativeubs[v], SCIPvarGetUbLocal(var), NULL, NULL) );
3362  (*branched) = TRUE;
3363
3364  SCIPdebugMsg(scip, "variable <%s> branched domain hole (%d,%g)\n", SCIPvarGetName(var),
3365  alternativeubs[v], SCIPvarGetUbLocal(var));
3366
3367  return SCIP_OKAY;
3368  }
3369  }
3370  }
3371
3372  return SCIP_OKAY;
3373 }
3374
3375 /** remove the capacity requirments for all job which start at the curtime */
3376 static
3378  SCIP_CONSDATA* consdata, /**< constraint data */
3379  int curtime, /**< current point in time */
3380  int* starttimes, /**< array of start times */
3381  int* startindices, /**< permutation with respect to the start times */
3382  int* freecapacity, /**< pointer to store the resulting free capacity */
3383  int* idx, /**< pointer to index in start time array */
3384  int nvars /**< number of vars in array of starttimes and startindices */
3385  )
3386 {
3387 #if defined SCIP_DEBUG && !defined NDEBUG
3388  int oldidx;
3389
3390  assert(idx != NULL);
3391  oldidx = *idx;
3392 #else
3393  assert(idx != NULL);
3394 #endif
3395
3396  assert(starttimes != NULL);
3397  assert(starttimes != NULL);
3398  assert(freecapacity != NULL);
3399  assert(starttimes[*idx] == curtime);
3400  assert(consdata->demands != NULL);
3401  assert(freecapacity != idx);
3402
3403  /* subtract all capacity needed up to this point */
3404  (*freecapacity) -= consdata->demands[startindices[*idx]];
3405
3406  while( (*idx)+1 < nvars && starttimes[(*idx)+1] == curtime )
3407  {
3408  ++(*idx);
3409  (*freecapacity) -= consdata->demands[startindices[(*idx)]];
3410  assert(freecapacity != idx);
3411  }
3412 #ifdef SCIP_DEBUG
3413  assert(oldidx <= *idx);
3414 #endif
3415 }
3416
3417 /** add the capacity requirments for all job which end at the curtime */
3418 static
3420  SCIP_CONSDATA* consdata, /**< constraint data */
3421  int curtime, /**< current point in time */
3422  int* endtimes, /**< array of end times */
3423  int* endindices, /**< permutation with rspect to the end times */
3424  int* freecapacity, /**< pointer to store the resulting free capacity */
3425  int* idx, /**< pointer to index in end time array */
3426  int nvars /**< number of vars in array of starttimes and startindices */
3427  )
3428 {
3429 #if defined SCIP_DEBUG && !defined NDEBUG
3430  int oldidx;
3431  oldidx = *idx;
3432 #endif
3433
3434  /* free all capacity usages of jobs the are no longer running */
3435  while( endtimes[*idx] <= curtime && *idx < nvars)
3436  {
3437  (*freecapacity) += consdata->demands[endindices[*idx]];
3438  ++(*idx);
3439  }
3440
3441 #ifdef SCIP_DEBUG
3442  assert(oldidx <= *idx);
3443 #endif
3444 }
3445
3446 /** computes a point in time when the capacity is exceeded returns hmax if this does not happen */
3447 static
3449  SCIP* scip, /**< SCIP data structure */
3450  SCIP_CONSDATA* consdata, /**< constraint handler data */
3451  SCIP_SOL* sol, /**< primal solution, or NULL for current LP/pseudo solution */
3452  int* timepoint /**< pointer to store the time point of the peak */
3453  )
3454 {
3455  int* starttimes; /* stores when each job is starting */
3456  int* endtimes; /* stores when each job ends */
3457  int* startindices; /* we will sort the startsolvalues, thus we need to know wich index of a job it corresponds to */
3458  int* endindices; /* we will sort the endsolvalues, thus we need to know wich index of a job it corresponds to */
3459
3460  int nvars; /* number of activities for this constraint */
3461  int freecapacity; /* remaining capacity */
3462  int curtime; /* point in time which we are just checking */
3463  int endindex; /* index of endsolvalues with: endsolvalues[endindex] > curtime */
3464
3465  int hmin;
3466  int hmax;
3467
3468  int j;
3469
3470  assert(consdata != NULL);
3471
3472  nvars = consdata->nvars;
3473  assert(nvars > 0);
3474
3475  *timepoint = consdata->hmax;
3476
3477  assert(consdata->vars != NULL);
3478
3479  SCIP_CALL( SCIPallocBufferArray(scip, &starttimes, nvars) );
3480  SCIP_CALL( SCIPallocBufferArray(scip, &endtimes, nvars) );
3481  SCIP_CALL( SCIPallocBufferArray(scip, &startindices, nvars) );
3482  SCIP_CALL( SCIPallocBufferArray(scip, &endindices, nvars) );
3483
3484  /* create event point arrays */
3485  createSortedEventpointsSol(scip, sol, consdata->nvars, consdata->vars, consdata->durations,
3486  starttimes, endtimes, startindices, endindices);
3487
3488  endindex = 0;
3489  freecapacity = consdata->capacity;
3490  hmin = consdata->hmin;
3491  hmax = consdata->hmax;
3492
3493  /* check each startpoint of a job whether the capacity is kept or not */
3494  for( j = 0; j < nvars; ++j )
3495  {
3496  curtime = starttimes[j];
3497  SCIPdebugMsg(scip, "look at %d-th job with start %d\n", j, curtime);
3498
3499  if( curtime >= hmax )
3500  break;
3501
3502  /* remove the capacity requirments for all job which start at the curtime */
3503  subtractStartingJobDemands(consdata, curtime, starttimes, startindices, &freecapacity, &j, nvars);
3504
3505  /* add the capacity requirments for all job which end at the curtime */
3506  addEndingJobDemands(consdata, curtime, endtimes, endindices, &freecapacity, &endindex, nvars);
3507
3508  assert(freecapacity <= consdata->capacity);
3509  assert(endindex <= nvars);
3510
3511  /* endindex - points to the next job which will finish */
3512  /* j - points to the last job that has been released */
3513
3514  /* if free capacity is smaller than zero, then add branching candidates */
3515  if( freecapacity < 0 && curtime >= hmin )
3516  {
3517  *timepoint = curtime;
3518  break;
3519  }
3520  } /*lint --e{850}*/
3521
3522  /* free all buffer arrays */
3523  SCIPfreeBufferArray(scip, &endindices);
3524  SCIPfreeBufferArray(scip, &startindices);
3525  SCIPfreeBufferArray(scip, &endtimes);
3526  SCIPfreeBufferArray(scip, &starttimes);
3527
3528  return SCIP_OKAY;
3529 }
3530
3531 /** checks all cumulative constraints for infeasibility and add branching candidates to storage */
3532 static
3534  SCIP* scip, /**< SCIP data structure */
3535  SCIP_CONS** conss, /**< constraints to be processed */
3536  int nconss, /**< number of constraints */
3537  SCIP_SOL* sol, /**< primal solution, or NULL for current LP/pseudo solution */
3538  int* nbranchcands /**< pointer to store the number of branching variables */
3539  )
3540 {
3541  SCIP_HASHTABLE* collectedvars;
3542  int c;
3543
3544  assert(scip != NULL);
3545  assert(conss != NULL);
3546
3547  /* create a hash table */
3548  SCIP_CALL( SCIPhashtableCreate(&collectedvars, SCIPblkmem(scip), SCIPgetNVars(scip),
3549  SCIPvarGetHashkey, SCIPvarIsHashkeyEq, SCIPvarGetHashkeyVal, NULL) );
3550
3551  assert(scip != NULL);
3552  assert(conss != NULL);
3553
3554  for( c = 0; c < nconss; ++c )
3555  {
3556  SCIP_CONS* cons;
3557  SCIP_CONSDATA* consdata;
3558
3559  int curtime;
3560  int j;
3561
3562  cons = conss[c];
3563  assert(cons != NULL);
3564
3565  if( !SCIPconsIsActive(cons) )
3566  continue;
3567
3568  consdata = SCIPconsGetData(cons);
3569  assert(consdata != NULL);
3570
3571  /* get point in time when capacity is exceeded */
3572  SCIP_CALL( computePeak(scip, consdata, sol, &curtime) );
3573
3574  if( curtime < consdata->hmin || curtime >= consdata->hmax )
3575  continue;
3576
3577  /* report all variables that are running at that point in time */
3578  for( j = 0; j < consdata->nvars; ++j )
3579  {
3580  SCIP_VAR* var;
3581  int lb;
3582  int ub;
3583
3584  var = consdata->vars[j];
3585  assert(var != NULL);
3586
3587  /* check if the variable was already added */
3588  if( SCIPhashtableExists(collectedvars, (void*)var) )
3589  continue;
3590
3591  lb = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
3592  ub = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
3593
3594  if( lb <= curtime && ub + consdata->durations[j] > curtime && lb < ub )
3595  {
3596  SCIP_Real solval;
3597  SCIP_Real score;
3598
3599  solval = SCIPgetSolVal(scip, sol, var);
3600  score = MIN(solval - lb, ub - solval) / ((SCIP_Real)ub-lb);
3601
3602  SCIPdebugMsg(scip, "add var <%s> to branch cand storage\n", SCIPvarGetName(var));
3603  SCIP_CALL( SCIPaddExternBranchCand(scip, var, score, lb + (ub - lb) / 2.0 + 0.2) );
3604  (*nbranchcands)++;
3605
3606  SCIP_CALL( SCIPhashtableInsert(collectedvars, var) );
3607  }
3608  }
3609  }
3610
3611  SCIPhashtableFree(&collectedvars);
3612
3613  SCIPdebugMsg(scip, "found %d branching candidates\n", *nbranchcands);
3614
3615  return SCIP_OKAY;
3616 }
3617
3618 /** enforcement of an LP, pseudo, or relaxation solution */
3619 static
3621  SCIP* scip, /**< SCIP data structure */
3622  SCIP_CONS** conss, /**< constraints to be processed */
3623  int nconss, /**< number of constraints */
3624  SCIP_SOL* sol, /**< solution to enforce (NULL for LP or pseudo solution) */
3625  SCIP_Bool branch, /**< should branching candidates be collected */
3626  SCIP_RESULT* result /**< pointer to store the result */
3627  )
3628 {
3629  if( branch )
3630  {
3631  int nbranchcands;
3632
3633  nbranchcands = 0;
3634  SCIP_CALL( collectBranchingCands(scip, conss, nconss, sol, &nbranchcands) );
3635
3636  if( nbranchcands > 0 )
3637  (*result) = SCIP_INFEASIBLE;
3638  }
3639  else
3640  {
3641  SCIP_Bool violated;
3642  int c;
3643
3644  violated = FALSE;
3645
3646  /* first check if a constraints is violated */
3647  for( c = 0; c < nconss && !violated; ++c )
3648  {
3649  SCIP_CONS* cons;
3650
3651  cons = conss[c];
3652  assert(cons != NULL);
3653
3654  SCIP_CALL( checkCons(scip, cons, sol, &violated, FALSE) );
3655  }
3656
3657  if( violated )
3658  (*result) = SCIP_INFEASIBLE;
3659  }
3660
3661  return SCIP_OKAY;
3662 }
3663
3664 /**@} */
3665
3666 /**@name Propagation
3667  *
3668  * @{
3669  */
3670
3671 /** check if cumulative constraint is independently of all other constraints */
3672 static
3674  SCIP* scip, /**< SCIP data structure */
3675  SCIP_CONS* cons /**< cumulative constraint */
3676  )
3677 {
3678  SCIP_CONSDATA* consdata;
3679  SCIP_VAR** vars;
3680  SCIP_Bool* downlocks;
3681  SCIP_Bool* uplocks;
3682  int nvars;
3683  int v;
3684
3685  consdata = SCIPconsGetData(cons);
3686  assert(consdata != NULL);
3687
3688  nvars = consdata->nvars;
3689  vars = consdata->vars;
3690  downlocks = consdata->downlocks;
3691  uplocks = consdata->uplocks;
3692
3693  /* check if the cumulative constraint has the only locks on the involved variables */
3694  for( v = 0; v < nvars; ++v )
3695  {
3696  SCIP_VAR* var;
3697
3698  var = vars[v];
3699  assert(var != NULL);
3700
3701  if( SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL) > (int)downlocks[v]
3702  || SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL) > (int)uplocks[v] )
3703  return FALSE;
3704  }
3705
3706  return TRUE;
3707 }
3708
3709 /** in case the cumulative constraint is independent of every else, solve the cumulative problem and apply the fixings
3710  * (dual reductions)
3711  */
3712 static
3714  SCIP* scip, /**< SCIP data structure */
3715  SCIP_CONS* cons, /**< cumulative constraint */
3716  SCIP_Longint maxnodes, /**< number of branch-and-bound nodes to solve an independent cumulative constraint (-1: no limit) */
3717  int* nchgbds, /**< pointer to store the number changed variable bounds */
3718  int* nfixedvars, /**< pointer to count number of fixings */
3719  int* ndelconss, /**< pointer to count number of deleted constraints */
3720  SCIP_Bool* cutoff, /**< pointer to store if the constraint is infeasible */
3721  SCIP_Bool* unbounded /**< pointer to store if the constraint is unbounded */
3722  )
3723 {
3724  SCIP_CONSDATA* consdata;
3725  SCIP_VAR** vars;
3726  SCIP_Real* objvals;
3727  SCIP_Real* lbs;
3728  SCIP_Real* ubs;
3729  SCIP_Real timelimit;
3730  SCIP_Real memorylimit;
3731  SCIP_Bool solved;
3732  SCIP_Bool error;
3733
3734  int ncheckconss;
3735  int nvars;
3736  int v;
3737
3738  assert(scip != NULL);
3739  assert(!SCIPconsIsModifiable(cons));
3740  assert(SCIPgetNConss(scip) > 0);
3741
3742  /* if SCIP is in probing mode or repropagation we cannot perform this dual reductions since this dual reduction
3743  * would/could end in an implication which can lead to cutoff of the/all optimal solution
3744  */
3745  if( SCIPinProbing(scip) || SCIPinRepropagation(scip) )
3746  return SCIP_OKAY;
3747
3748  /* constraints for which the check flag is set to FALSE, did not contribute to the lock numbers; therefore, we cannot
3749  * use the locks to decide for a dual reduction using this constraint;
3750  */
3751  if( !SCIPconsIsChecked(cons) )
3752  return SCIP_OKAY;
3753
3754  ncheckconss = SCIPgetNCheckConss(scip);
3755
3756  /* if the cumulative constraint is the only constraint of the original problem or the only check constraint in the
3757  * presolved problem do nothing execpt to change the parameter settings
3758  */
3759  if( ncheckconss == 1 )
3760  {
3761  /* shrink the minimal maximum value for the conflict length */
3762  SCIP_CALL( SCIPsetIntParam(scip, "conflict/minmaxvars", 10) );
3763
3764  /* use only first unique implication point */
3765  SCIP_CALL( SCIPsetIntParam(scip, "conflict/fuiplevels", 1) );
3766
3767  /* do not use reconversion conflicts */
3768  SCIP_CALL( SCIPsetIntParam(scip, "conflict/reconvlevels", 0) );
3769
3770  /* after 250 conflict we force a restart since then the variable statistics are reasonable initialized */
3771  SCIP_CALL( SCIPsetIntParam(scip, "conflict/restartnum", 250) );
3772
3773  /* increase the number of conflicts which induce a restart */
3774  SCIP_CALL( SCIPsetRealParam(scip, "conflict/restartfac", 2.0) );
3775
3776  /* weight the variable which made into a conflict */
3777  SCIP_CALL( SCIPsetRealParam(scip, "conflict/conflictweight", 1.0) );
3778
3779  /* do not check pseudo solution (for performance reasons) */
3780  SCIP_CALL( SCIPsetBoolParam(scip, "constraints/disableenfops", TRUE) );
3781
3782  /* use value based history to detect a reasonable branching point */
3783  SCIP_CALL( SCIPsetBoolParam(scip, "history/valuebased", TRUE) );
3784
3785  /* turn of LP relaxation */
3786  SCIP_CALL( SCIPsetIntParam(scip, "lp/solvefreq", -1) );
3787
3788  /* prefer the down branch in case the value based history does not suggest something */
3789  SCIP_CALL( SCIPsetCharParam(scip, "nodeselection/childsel", 'd') );
3790
3791  /* accept any bound change */
3792  SCIP_CALL( SCIPsetRealParam(scip, "numerics/boundstreps", 1e-6) );
3793
3794  /* allow for at most 10 restart, after that the value based history should be reliable */
3795  SCIP_CALL( SCIPsetIntParam(scip, "presolving/maxrestarts", 10) );
3796
3797  /* set priority for depth first search to highest possible value */
3798  SCIP_CALL( SCIPsetIntParam(scip, "nodeselection/dfs/stdpriority", INT_MAX/4) );
3799
3800  return SCIP_OKAY;
3801  }
3802
3803  consdata = SCIPconsGetData(cons);
3804  assert(consdata != NULL);
3805
3806  /* check if already tried to solve that constraint as independent sub problem; we do not want to try it again if we
3807  * fail on the first place
3808  */
3809  if( consdata->triedsolving )
3810  return SCIP_OKAY;
3811
3812  /* check if constraint is independently */
3813  if( !isConsIndependently(scip, cons) )
3814  return SCIP_OKAY;
3815
3816  /* mark the constraint to be tried of solving it as independent sub problem; in case that is successful the
3817  * constraint is deleted; otherwise, we want to ensure that we do not try that again
3818  */
3819  consdata->triedsolving = TRUE;
3820
3821  SCIPdebugMsg(scip, "the cumulative constraint <%s> is independent from rest of the problem (%d variables, %d constraints)\n",
3822  SCIPconsGetName(cons), SCIPgetNVars(scip), SCIPgetNConss(scip));
3823  SCIPdebugPrintCons(scip, cons, NULL);
3824
3825  nvars = consdata->nvars;
3826  vars = consdata->vars;
3827
3828  SCIP_CALL( SCIPallocBufferArray(scip, &lbs, nvars) );
3829  SCIP_CALL( SCIPallocBufferArray(scip, &ubs, nvars) );
3830  SCIP_CALL( SCIPallocBufferArray(scip, &objvals, nvars) );
3831
3832  for( v = 0; v < nvars; ++v )
3833  {
3834  SCIP_VAR* var;
3835
3836  /* if a variables array is given, use the variable bounds otherwise the default values stored in the ests and lsts
3837  * array
3838  */
3839  var = vars[v];
3840  assert(var != NULL);
3841
3842  lbs[v] = SCIPvarGetLbLocal(var);
3843  ubs[v] = SCIPvarGetUbLocal(var);
3844
3845  objvals[v] = SCIPvarGetObj(var);
3846  }
3847
3848  /* check whether there is enough time and memory left */
3849  SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
3850  if( !SCIPisInfinity(scip, timelimit) )
3851  timelimit -= SCIPgetSolvingTime(scip);
3852  SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &memorylimit) );
3853
3854  /* substract the memory already used by the main SCIP and the estimated memory usage of external software */
3855  if( !SCIPisInfinity(scip, memorylimit) )
3856  {
3857  memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
3858  memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
3859  }
3860
3861  /* solve the cumulative condition separately */
3862  SCIP_CALL( SCIPsolveCumulative(scip, nvars, lbs, ubs, objvals, consdata->durations, consdata->demands, consdata->capacity,
3863  consdata->hmin, consdata->hmax, timelimit, memorylimit, maxnodes, &solved, cutoff, unbounded, &error) );
3864
3865  if( !(*cutoff) && !(*unbounded) && !error )
3866  {
3867  SCIP_Bool infeasible;
3868  SCIP_Bool tightened;
3869  SCIP_Bool allfixed;
3870
3871  allfixed = TRUE;
3872
3873  for( v = 0; v < nvars; ++v )
3874  {
3875  /* check if variable is fixed */
3876  if( lbs[v] + 0.5 > ubs[v] )
3877  {
3878  SCIP_CALL( SCIPfixVar(scip, vars[v], lbs[v], &infeasible, &tightened) );
3879  assert(!infeasible);
3880
3881  if( tightened )
3882  {
3883  (*nfixedvars)++;
3884  consdata->triedsolving = FALSE;
3885  }
3886  }
3887  else
3888  {
3889  SCIP_CALL( SCIPtightenVarLb(scip, vars[v], lbs[v], TRUE, &infeasible, &tightened) );
3890  assert(!infeasible);
3891
3892  if( tightened )
3893  {
3894  (*nchgbds)++;
3895  consdata->triedsolving = FALSE;
3896  }
3897
3898  SCIP_CALL( SCIPtightenVarUb(scip, vars[v], ubs[v], TRUE, &infeasible, &tightened) );
3899  assert(!infeasible);
3900
3901  if( tightened )
3902  {
3903  (*nchgbds)++;
3904  consdata->triedsolving = FALSE;
3905  }
3906
3907  allfixed = FALSE;
3908  }
3909  }
3910
3911  /* if all variables are fixed, remove the cumulative constraint since it is redundant */
3912  if( allfixed )
3913  {
3914  SCIP_CALL( SCIPdelConsLocal(scip, cons) );
3915  (*ndelconss)++;
3916  }
3917  }
3918
3919  SCIPfreeBufferArray(scip, &objvals);
3920  SCIPfreeBufferArray(scip, &ubs);
3921  SCIPfreeBufferArray(scip, &lbs);
3922
3923  return SCIP_OKAY;
3924 }
3925
3926 /** start conflict analysis to analysis the core insertion which is infeasible */
3927 static
3929  SCIP* scip, /**< SCIP data structure */
3930  int nvars, /**< number of start time variables (activities) */
3931  SCIP_VAR** vars, /**< array of start time variables */
3932  int* durations, /**< array of durations */
3933  int* demands, /**< array of demands */
3934  int capacity, /**< cumulative capacity */
3935  int hmin, /**< left bound of time axis to be considered (including hmin) */
3936  int hmax, /**< right bound of time axis to be considered (not including hmax) */
3937  SCIP_VAR* infervar, /**< start time variable which lead to the infeasibilty */
3938  int inferduration, /**< duration of the start time variable */
3939  int inferdemand, /**< demand of the start time variable */
3940  int inferpeak, /**< profile preak which causes the infeasibilty */
3941  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
3942  SCIP_Bool* initialized, /**< pointer to store if the conflict analysis was initialized */
3943  SCIP_Bool* explanation /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
3944  )
3945 {
3946  SCIPdebugMsg(scip, "detected infeasibility due to adding a core to the core resource profile\n");
3947  SCIPdebugMsg(scip, "variable <%s>[%g,%g] (demand %d, duration %d)\n", SCIPvarGetName(infervar),
3948  SCIPvarGetLbLocal(infervar), SCIPvarGetUbLocal(infervar), inferdemand, inferduration);
3949
3950  /* initialize conflict analysis if conflict analysis is applicable */
3952  {
3954
3955  SCIP_CALL( resolvePropagationCoretimes(scip, nvars, vars, durations, demands, capacity, hmin, hmax,
3956  infervar, inferdemand, inferpeak, inferpeak, NULL, usebdwidening, NULL, explanation) );
3957
3958  SCIPdebugMsg(scip, "add lower and upper bounds of variable <%s>\n", SCIPvarGetName(infervar));
3959
3960  /* add both bound of the inference variable since these biuld the core which we could not inserted */
3961  if( usebdwidening )
3962  {
3963  SCIP_CALL( SCIPaddConflictRelaxedLb(scip, infervar, NULL, (SCIP_Real)(inferpeak - inferduration + 1)) );
3964  SCIP_CALL( SCIPaddConflictRelaxedUb(scip, infervar, NULL, (SCIP_Real)inferpeak) );
3965  }
3966  else
3967  {
3968  SCIP_CALL( SCIPaddConflictLb(scip, infervar, NULL) );
3969  SCIP_CALL( SCIPaddConflictUb(scip, infervar, NULL) );
3970  }
3971
3972  *initialized = TRUE;
3973  }
3974
3975  return SCIP_OKAY;
3976 }
3977
3978 /** We are using the core resource profile which contains all core except the one of the start time variable which we
3979  * want to propagate, to incease the earliest start time. This we are doing in steps of length at most the duration of
3980  * the job. The reason for that is, that this makes it later easier to resolve this propagation during the conflict
3981  * analysis
3982  */
3983 static
3985  SCIP* scip, /**< SCIP data structure */
3986  int nvars, /**< number of start time variables (activities) */
3987  SCIP_VAR** vars, /**< array of start time variables */
3988  int* durations, /**< array of durations */
3989  int* demands, /**< array of demands */
3990  int capacity, /**< cumulative capacity */
3991  int hmin, /**< left bound of time axis to be considered (including hmin) */
3992  int hmax, /**< right bound of time axis to be considered (not including hmax) */
3993  SCIP_CONS* cons, /**< constraint which is propagated */
3994  SCIP_PROFILE* profile, /**< resource profile */
3995  int idx, /**< position of the variable to propagate */
3996  int* nchgbds, /**< pointer to store the number of bound changes */
3997  SCIP_Bool usebdwidening, /**< should bound widening be used during conflict analysis? */
3998  SCIP_Bool* initialized, /**< was conflict analysis initialized */
3999  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
4000  SCIP_Bool* infeasible /**< pointer to store if the constraint is infeasible */
4001  )
4002 {
4003  SCIP_VAR* var;
4004  int ntimepoints;
4005  int duration;
4006  int demand;
4007  int peak;
4008  int newlb;
4009  int est;
4010  int lst;
4011  int pos;
4012
4013  var = vars[idx];
4014  assert(var != NULL);
4015
4016  duration = durations[idx];
4017  assert(duration > 0);
4018
4019  demand = demands[idx];
4020  assert(demand > 0);
4021
4022  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
4023  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
4024  ntimepoints = SCIPprofileGetNTimepoints(profile);
4025
4026  /* first we find left position of earliest start time (lower bound) in resource profile; this position gives us the
4027  * load which we have at the earliest start time (lower bound)
4028  */
4029  (void) SCIPprofileFindLeft(profile, est, &pos);
4030
4031  SCIPdebugMsg(scip, "propagate earliest start time (lower bound) (pos %d)\n", pos);
4032
4033  /* we now trying to move the earliest start time in steps of at most "duration" length */
4034  do
4035  {
4036  INFERINFO inferinfo;
4037  SCIP_Bool tightened;
4038  int ect;
4039
4040 #ifndef NDEBUG
4041  {
4042  /* in debug mode we check that we adjust the search position correctly */
4043  int tmppos;
4044
4045  (void)SCIPprofileFindLeft(profile, est, &tmppos);
4046  assert(pos == tmppos);
4047  }
4048 #endif
4049  ect = est + duration;
4050  peak = -1;
4051
4052  /* we search for a peak within the core profile which conflicts with the demand of the start time variable; we
4053  * want a peak which is closest to the earliest completion time
4054  */
4055  do
4056  {
4057  /* check if the profile load conflicts with the demand of the start time variable */
4058  if( SCIPprofileGetLoad(profile, pos) + demand > capacity )
4059  peak = pos;
4060
4061  pos++;
4062  }
4063  while( pos < ntimepoints && SCIPprofileGetTime(profile, pos) < ect );
4064
4065  /* if we found no peak that means current the job could be scheduled at its earliest start time without
4066  * conflicting to the core resource profile
4067  */
4068  if( peak == -1 )
4069  break;
4070
4071  /* the peak position gives us a time point where the start time variable is in conflict with the resource
4072  * profile. That means we have to move it to the next time point in the resource profile but at most to the
4073  * earliest completion time (the remaining move will done in the next loop)
4074  */
4075  newlb = SCIPprofileGetTime(profile, peak+1);
4076  newlb = MIN(newlb, ect);
4077
4078  /* if the earliest start time is greater than the lst we detected an infeasibilty */
4079  if( newlb > lst )
4080  {
4081  SCIPdebugMsg(scip, "variable <%s>: cannot be scheduled\n", SCIPvarGetName(var));
4082
4083  /* use conflict analysis to analysis the core insertion which was infeasible */
4084  SCIP_CALL( analyseInfeasibelCoreInsertion(scip, nvars, vars, durations, demands, capacity, hmin, hmax,
4085  var, duration, demand, newlb-1, usebdwidening, initialized, explanation) );
4086
4087  if( explanation != NULL )
4088  explanation[idx] = TRUE;
4089
4090  *infeasible = TRUE;
4091
4092  break;
4093  }
4094
4095  /* construct the inference information which we are using with the conflict analysis to resolve that particular
4096  * bound change
4097  */
4098  inferinfo = getInferInfo(PROPRULE_1_CORETIMES, idx, newlb-1);
4099
4100  /* perform the bound lower bound change */
4101  SCIP_CALL( SCIPinferVarLbCons(scip, var, (SCIP_Real)newlb, cons, inferInfoToInt(inferinfo), TRUE, infeasible, &tightened) );
4102  assert(tightened);
4103  assert(!(*infeasible));
4104
4105  SCIPdebugMsg(scip, "variable <%s> new lower bound <%d> -> <%d>\n", SCIPvarGetName(var), est, newlb);
4106  (*nchgbds)++;
4107
4108  /* for the statistic we count the number of times a lower bound was tightened due the the time-table algorithm */
4110
4111  /* adjust the earliest start time
4112  *
4113  * @note We are taking the lower of the start time variable on purpose instead of newlb. This is due the fact that
4114  * the proposed lower bound might be even strength by be the core which can be the case if aggregations are
4115  * involved.
4116  */
4117  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
4118  assert(est >= newlb);
4119
4120  /* adjust the search position for the resource profile for the next step */
4121  if( est == SCIPprofileGetTime(profile, peak+1) )
4122  pos = peak + 1;
4123  else
4124  pos = peak;
4125  }
4126  while( est < lst );
4127
4128  return SCIP_OKAY;
4129 }
4130
4131 /** We are using the core resource profile which contains all core except the one of the start time variable which we
4132  * want to propagate, to decrease the latest start time. This we are doing in steps of length at most the duration of
4133  * the job. The reason for that is, that this makes it later easier to resolve this propagation during the conflict
4134  * analysis
4135  */
4136 static
4138  SCIP* scip, /**< SCIP data structure */
4139  SCIP_VAR* var, /**< start time variable to propagate */
4140  int duration, /**< duration of the job */
4141  int demand, /**< demand of the job */
4142  int capacity, /**< cumulative capacity */
4143  SCIP_CONS* cons, /**< constraint which is propagated */
4144  SCIP_PROFILE* profile, /**< resource profile */
4145  int idx, /**< position of the variable to propagate */
4146  int* nchgbds /**< pointer to store the number of bound changes */
4147  )
4148 {
4149  int ntimepoints;
4150  int newub;
4151  int peak;
4152  int pos;
4153  int est;
4154  int lst;
4155  int lct;
4156
4157  assert(var != NULL);
4158  assert(duration > 0);
4159  assert(demand > 0);
4160
4161  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
4162  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
4163
4164  /* in case the start time variable is fixed do nothing */
4165  if( est == lst )
4166  return SCIP_OKAY;
4167
4168  ntimepoints = SCIPprofileGetNTimepoints(profile);
4169
4170  lct = lst + duration;
4171
4172  /* first we find left position of latest completion time minus 1 (upper bound + duration) in resource profile; That
4173  * is the last time point where the job would run if schedule it at its latest start time (upper bound). This
4174  * position gives us the load which we have at the latest completion time minus one
4175  */
4176  (void) SCIPprofileFindLeft(profile, lct - 1, &pos);
4177
4178  SCIPdebugMsg(scip, "propagate upper bound (pos %d)\n", pos);
4179  SCIPdebug( SCIPprofilePrint(profile, SCIPgetMessagehdlr(scip), NULL) );
4180
4181  if( pos == ntimepoints-1 && SCIPprofileGetTime(profile, pos) == lst )
4182  return SCIP_OKAY;
4183
4184  /* we now trying to move the latest start time in steps of at most "duration" length */
4185  do
4186  {
4187  INFERINFO inferinfo;
4188  SCIP_Bool tightened;
4189  SCIP_Bool infeasible;
4190
4191  peak = -1;
4192
4193 #ifndef NDEBUG
4194  {
4195  /* in debug mode we check that we adjust the search position correctly */
4196  int tmppos;
4197
4198  (void)SCIPprofileFindLeft(profile, lct - 1, &tmppos);
4199  assert(pos == tmppos);
4200  }
4201 #endif
4202
4203  /* we search for a peak within the core profile which conflicts with the demand of the start time variable; we
4204  * want a peak which is closest to the latest start time
4205  */
4206  do
4207  {
4208  if( SCIPprofileGetLoad(profile, pos) + demand > capacity )
4209  peak = pos;
4210
4211  pos--;
4212  }
4213  while( pos >= 0 && SCIPprofileGetTime(profile, pos+1) > lst);
4214
4215  /* if we found no peak that means the current job could be scheduled at its latest start time without conflicting
4216  * to the core resource profile
4217  */
4218  if( peak == -1 )
4219  break;
4220
4221  /* the peak position gives us a time point where the start time variable is in conflict with the resource
4222  * profile. That means the job has be done until that point. Hence that gives us the latest completion
4223  * time. Note that that we want to move the bound by at most the duration length (the remaining move we are
4224  * doing in the next loop)
4225  */
4226  newub = SCIPprofileGetTime(profile, peak);
4227  newub = MAX(newub, lst) - duration;
4228  assert(newub >= est);
4229
4230  /* construct the inference information which we are using with the conflict analysis to resolve that particular
4231  * bound change
4232  */
4233  inferinfo = getInferInfo(PROPRULE_1_CORETIMES, idx, newub+duration);
4234
4235  /* perform the bound upper bound change */
4236  SCIP_CALL( SCIPinferVarUbCons(scip, var, (SCIP_Real)newub, cons, inferInfoToInt(inferinfo), TRUE, &infeasible, &tightened) );
4237  assert(tightened);
4238  assert(!infeasible);
4239
4240  SCIPdebugMsg(scip, "variable <%s>: new upper bound <%d> -> <%d>\n", SCIPvarGetName(var), lst, newub);
4241  (*nchgbds)++;
4242
4243  /* for the statistic we count the number of times a upper bound was tightened due the the time-table algorithm */
4245
4246  /* adjust the latest start and completion time
4247  *
4248  * @note We are taking the upper of the start time variable on purpose instead of newub. This is due the fact that
4249  * the proposed upper bound might be even strength by be the core which can be the case if aggregations are
4250  * involved.
4251  */
4252  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
4253  assert(lst <= newub);
4254  lct = lst + duration;
4255
4256  /* adjust the search position for the resource profile for the next step */
4257  if( SCIPprofileGetTime(profile, peak) == lct )
4258  pos = peak - 1;
4259  else
4260  pos = peak;
4261  }
4262  while( est < lst );
4263
4264  return SCIP_OKAY;
4265 }
4266
4267 /** compute for the different earliest start and latest completion time the core energy of the corresponding time
4268  * points
4269  */
4270 static
4272  SCIP* scip, /**< SCIP data structure */
4273  SCIP_PROFILE* profile, /**< core profile */
4274  int nvars, /**< number of start time variables (activities) */
4275  int* ests, /**< array of sorted earliest start times */
4276  int* lcts, /**< array of sorted latest completion times */
4277  int* coreEnergyAfterEst, /**< array to store the core energy after the earliest start time of each job */
4278  int* coreEnergyAfterLct /**< array to store the core energy after the latest completion time of each job */
4279  )
4280 {
4281  int ntimepoints;
4282  int energy;
4283  int t;
4284  int v;
4285
4286  ntimepoints = SCIPprofileGetNTimepoints(profile);
4287  t = ntimepoints - 1;
4288  energy = 0;
4289
4290  /* compute core energy after the earliest start time of each job */
4291  for( v = nvars-1; v >= 0; --v )
4292  {
4293  while( t > 0 && SCIPprofileGetTime(profile, t-1) >= ests[v] )
4294  {
4295  assert(SCIPprofileGetLoad(profile, t-1) >= 0);
4296  assert(SCIPprofileGetTime(profile, t) - SCIPprofileGetTime(profile, t-1)>= 0);
4297  energy += SCIPprofileGetLoad(profile, t-1) * (SCIPprofileGetTime(profile, t) - SCIPprofileGetTime(profile, t-1));
4298  t--;
4299  }
4300  assert(SCIPprofileGetTime(profile, t) >= ests[v] || t == ntimepoints-1);
4301
4302  /* maybe ests[j] is in-between two timepoints */
4303  if( SCIPprofileGetTime(profile, t) - ests[v] > 0 )
4304  {
4305  assert(t > 0);
4306  coreEnergyAfterEst[v] = energy + SCIPprofileGetLoad(profile, t-1) * (SCIPprofileGetTime(profile, t) - ests[v]);
4307  }
4308  else
4309  coreEnergyAfterEst[v] = energy;
4310  }
4311
4312  t = ntimepoints - 1;
4313  energy = 0;
4314
4315  /* compute core energy after the latest completion time of each job */
4316  for( v = nvars-1; v >= 0; --v )
4317  {
4318  while( t > 0 && SCIPprofileGetTime(profile, t-1) >= lcts[v] )
4319  {
4320  assert(SCIPprofileGetLoad(profile, t-1) >= 0);
4321  assert(SCIPprofileGetTime(profile, t) - SCIPprofileGetTime(profile, t-1)>= 0);
4322  energy += SCIPprofileGetLoad(profile, t-1) * (SCIPprofileGetTime(profile, t) - SCIPprofileGetTime(profile, t-1));
4323  t--;
4324  }
4325  assert(SCIPprofileGetTime(profile, t) >= lcts[v] || t == ntimepoints-1);
4326
4327  /* maybe lcts[j] is in-between two timepoints */
4328  if( SCIPprofileGetTime(profile, t) - lcts[v] > 0 )
4329  {
4330  assert(t > 0);
4331  coreEnergyAfterLct[v] = energy + SCIPprofileGetLoad(profile, t-1) * (SCIPprofileGetTime(profile, t) - lcts[v]);
4332  }
4333  else
4334  coreEnergyAfterLct[v] = energy;
4335  }
4336
4337  return SCIP_OKAY;
4338 }
4339
4340 /** collect earliest start times, latest completion time, and free energy contributions */
4341 static
4342 void collectDataTTEF(
4343  SCIP* scip, /**< SCIP data structure */
4344  int nvars, /**< number of start time variables (activities) */
4345  SCIP_VAR** vars, /**< array of start time variables */
4346  int* durations, /**< array of durations */
4347  int* demands, /**< array of demands */
4348  int hmin, /**< left bound of time axis to be considered (including hmin) */
4349  int hmax, /**< right bound of time axis to be considered (not including hmax) */
4350  int* permests, /**< array to store the variable positions */
4351  int* ests, /**< array to store earliest start times */
4352  int* permlcts, /**< array to store the variable positions */
4353  int* lcts, /**< array to store latest completion times */
4354  int* ects, /**< array to store earliest completion times of the flexible part of the job */
4355  int* lsts, /**< array to store latest start times of the flexible part of the job */
4356  int* flexenergies /**< array to store the flexible energies of each job */
4357  )
4358 {
4359  int v;
4360
4361  for( v = 0; v < nvars; ++ v)
4362  {
4363  int duration;
4366  int core;
4367  int est;
4368  int lct;
4369  int ect;
4370  int lst;
4371
4372  duration = durations[v];
4373  assert(duration > 0);
4374
4375  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[v]));
4376  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[v]));
4377  ect = est + duration;
4378  lct = lst + duration;
4379
4380  ests[v] = est;
4381  lcts[v] = lct;
4382  permests[v] = v;
4383  permlcts[v] = v;
4384
4385  /* compute core time window which lies within the effective horizon */
4386  core = (int) computeCoreWithInterval(hmin, hmax, ect, lst);
4387
4388  /* compute the number of time steps the job could run before the effective horizon */
4389  leftadjust = MAX(0, hmin - est);
4390
4391  /* compute the number of time steps the job could run after the effective horizon */
4392  rightadjust = MAX(0, lct - hmax);
4393
4394  /* compute for each job the energy which is flexible; meaning not part of the core */
4395  flexenergies[v] = duration - leftadjust - rightadjust - core;
4396  flexenergies[v] = MAX(0, flexenergies[v]);
4397  flexenergies[v] *= demands[v];
4398  assert(flexenergies[v] >= 0);
4399
4400  /* the earliest completion time of the flexible energy */
4401  ects[v] = MIN(ect, lst);
4402
4403  /* the latest start time of the flexible energy */
4404  lsts[v] = MAX(ect, lst);
4405  }
4406 }
4407
4408 /** try to tighten the lower bound of the given variable */
4409 static
4411  SCIP* scip, /**< SCIP data structure */
4412  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
4413  int nvars, /**< number of start time variables (activities) */
4414  SCIP_VAR** vars, /**< array of start time variables */
4415  int* durations, /**< array of durations */
4416  int* demands, /**< array of demands */
4417  int capacity, /**< cumulative capacity */
4418  int hmin, /**< left bound of time axis to be considered (including hmin) */
4419  int hmax, /**< right bound of time axis to be considered (not including hmax) */
4420  SCIP_VAR* var, /**< variable to be considered for upper bound tightening */
4421  int duration, /**< duration of the job */
4422  int demand, /**< demand of the job */
4423  int est, /**< earliest start time of the job */
4424  int ect, /**< earliest completion time of the flexible part of the job */
4425  int lct, /**< latest completion time of the job */
4426  int begin, /**< begin of the time window under investigation */
4427  int end, /**< end of the time window under investigation */
4428  SCIP_Longint energy, /**< available energy for the flexible part of the hob within the time window */
4429  int* bestlb, /**< pointer to strope the best lower bound change */
4430  int* inferinfos, /**< pointer to store the inference information which is need for the (best) lower bound change */
4431  SCIP_Bool* initialized, /**< was conflict analysis initialized */
4432  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
4433  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
4434  )
4435 {
4436  int newlb;
4437
4438  assert(begin >= hmin);
4439  assert(end <= hmax);
4440
4441  /* check if the time-table edge-finding should infer bounds */
4442  if( !conshdlrdata->ttefinfer )
4443  return SCIP_OKAY;
4444
4445  /* if the job can be processed completely before or after the time window, nothing can be tightened */
4446  if( est >= end || ect <= begin )
4447  return SCIP_OKAY;
4448
4449  /* if flexible part runs completely within the time window (assuming it is scheduled on its earliest start time), we
4450  * skip since the overload check will do the job
4451  */
4452  if( est >= begin && ect <= end )
4453  return SCIP_OKAY;
4454
4455  /* check if the available energy in the time window is to small to handle the flexible part if it is schedule on its
4456  * earliest start time
4457  */
4458  if( energy >= demand * ((SCIP_Longint) MAX(begin, est) - MIN(end, ect)) )
4459  return SCIP_OKAY;
4460
4461  /* adjust the available energy for the job; the given available energy assumes that the core of the considered job is
4462  * present; therefore, we need to add the core;
4463  *
4464  * @note the variable ect define the earliest completion time of the flexible part of the job; hence we need to
4465  * compute the earliest completion time of the (whole) job
4466  */
4467  energy += computeCoreWithInterval(begin, end, est + duration, lct - duration) * demand;
4468
4469  /* compute a latest start time (upper bound) such that the job consums at most the available energy
4470  *
4471  * @note we can round down the compute duration w.r.t. the available energy
4472  */
4473  newlb = end - (int) (energy / demand);
4474
4475  /* check if we detected an infeasibility which is the case if the new lower bound is larger than the current upper
4476  * bound (latest start time); meaning it is not possible to schedule the job
4477  */
4478  if( newlb > lct - duration )
4479  {
4480  /* initialize conflict analysis if conflict analysis is applicable */
4482  {
4483  SCIP_Real relaxedbd;
4484
4485  assert(SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var)) < newlb);
4486
4487  /* it is enough to overshoot the upper bound of the variable by one */
4488  relaxedbd = SCIPvarGetUbLocal(var) + 1.0;
4489
4490  /* initialize conflict analysis */
4492
4493  /* added to upper bound (which was overcut be new lower bound) of the variable */
4494  SCIP_CALL( SCIPaddConflictUb(scip, var, NULL) );
4495
4496  /* analyze the infeasible */
4497  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
4498  begin, end, var, SCIP_BOUNDTYPE_LOWER, NULL, relaxedbd, conshdlrdata->usebdwidening, explanation) );
4499
4500  (*initialized) = TRUE;
4501  }
4502
4503  (*cutoff) = TRUE;
4504  }
4505  else if( newlb > (*bestlb) )
4506  {
4507  INFERINFO inferinfo;
4508
4509  assert(newlb > begin);
4510
4511  inferinfo = getInferInfo(PROPRULE_3_TTEF, begin, end);
4512
4513  /* construct inference information */
4514  (*inferinfos) = inferInfoToInt(inferinfo);
4515  (*bestlb) = newlb;
4516  }
4517
4518  return SCIP_OKAY;
4519 }
4520
4521 /** try to tighten the upper bound of the given variable */
4522 static
4524  SCIP* scip, /**< SCIP data structure */
4525  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
4526  int nvars, /**< number of start time variables (activities) */
4527  SCIP_VAR** vars, /**< array of start time variables */
4528  int* durations, /**< array of durations */
4529  int* demands, /**< array of demands */
4530  int capacity, /**< cumulative capacity */
4531  int hmin, /**< left bound of time axis to be considered (including hmin) */
4532  int hmax, /**< right bound of time axis to be considered (not including hmax) */
4533  SCIP_VAR* var, /**< variable to be considered for upper bound tightening */
4534  int duration, /**< duration of the job */
4535  int demand, /**< demand of the job */
4536  int est, /**< earliest start time of the job */
4537  int lst, /**< latest start time of the flexible part of the job */
4538  int lct, /**< latest completion time of the job */
4539  int begin, /**< begin of the time window under investigation */
4540  int end, /**< end of the time window under investigation */
4541  SCIP_Longint energy, /**< available energy for the flexible part of the hob within the time window */
4542  int* bestub, /**< pointer to strope the best upper bound change */
4543  int* inferinfos, /**< pointer to store the inference information which is need for the (best) upper bound change */
4544  SCIP_Bool* initialized, /**< was conflict analysis initialized */
4545  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
4546  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
4547  )
4548 {
4549  int newub;
4550
4551  assert(begin >= hmin);
4552  assert(end <= hmax);
4553  assert(est < begin);
4554
4555  /* check if the time-table edge-finding should infer bounds */
4556  if( !conshdlrdata->ttefinfer )
4557  return SCIP_OKAY;
4558
4559  /* if flexible part of the job can be processed completely before or after the time window, nothing can be tightened */
4560  if( lst >= end || lct <= begin )
4561  return SCIP_OKAY;
4562
4563  /* if flexible part runs completely within the time window (assuming it is scheduled on its latest start time), we
4564  * skip since the overload check will do the job
4565  */
4566  if( lst >= begin && lct <= end )
4567  return SCIP_OKAY;
4568
4569  /* check if the available energy in the time window is to small to handle the flexible part of the job */
4570  if( energy >= demand * ((SCIP_Longint) MIN(end, lct) - MAX(begin, lst)) )
4571  return SCIP_OKAY;
4572
4573  /* adjust the available energy for the job; the given available energy assumes that the core of the considered job is
4574  * present; therefore, we need to add the core;
4575  *
4576  * @note the variable lst define the latest start time of the flexible part of the job; hence we need to compute the
4577  * latest start of the (whole) job
4578  */
4579  energy += computeCoreWithInterval(begin, end, est + duration, lct - duration) * demand;
4580  assert(energy >= 0);
4581
4582  /* compute a latest start time (upper bound) such that the job consums at most the available energy
4583  *
4584  * @note we can round down the compute duration w.r.t. the available energy
4585  */
4586  assert(demand > 0);
4587  newub = begin - duration + (int) (energy / demand);
4588
4589  /* check if we detected an infeasibility which is the case if the new upper bound is smaller than the current lower
4590  * bound (earliest start time); meaning it is not possible to schedule the job
4591  */
4592  if( newub < est )
4593  {
4594  /* initialize conflict analysis if conflict analysis is applicable */
4596  {
4597  SCIP_Real relaxedbd;
4598
4599  assert(SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var)) > newub);
4600
4601  /* it is enough to undershoot the lower bound of the variable by one */
4602  relaxedbd = SCIPvarGetLbLocal(var) - 1.0;
4603
4604  /* initialize conflict analysis */
4606
4607  /* added to lower bound (which was undercut be new upper bound) of the variable */
4608  SCIP_CALL( SCIPaddConflictLb(scip, var, NULL) );
4609
4610  /* analyze the infeasible */
4611  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
4612  begin, end, var, SCIP_BOUNDTYPE_UPPER, NULL, relaxedbd, conshdlrdata->usebdwidening, explanation) );
4613
4614  (*initialized) = TRUE;
4615  }
4616
4617  (*cutoff) = TRUE;
4618  }
4619  else if( newub < (*bestub) )
4620  {
4621  INFERINFO inferinfo;
4622
4623  assert(newub < begin);
4624
4625  inferinfo = getInferInfo(PROPRULE_3_TTEF, begin, end);
4626
4627  /* construct inference information */
4628  (*inferinfos) = inferInfoToInt(inferinfo);
4629  (*bestub) = newub;
4630  }
4631
4632  return SCIP_OKAY;
4633 }
4634
4635 /** propagate the upper bounds and "opportunistically" the lower bounds using the time-table edge-finding algorithm */
4636 static
4638  SCIP* scip, /**< SCIP data structure */
4639  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
4640  int nvars, /**< number of start time variables (activities) */
4641  SCIP_VAR** vars, /**< array of start time variables */
4642  int* durations, /**< array of durations */
4643  int* demands, /**< array of demands */
4644  int capacity, /**< cumulative capacity */
4645  int hmin, /**< left bound of time axis to be considered (including hmin) */
4646  int hmax, /**< right bound of time axis to be considered (not including hmax) */
4647  int* newlbs, /**< array to buffer new lower bounds */
4648  int* newubs, /**< array to buffer new upper bounds */
4649  int* lbinferinfos, /**< array to store the inference information for the lower bound changes */
4650  int* ubinferinfos, /**< array to store the inference information for the upper bound changes */
4651  int* lsts, /**< array of latest start time of the flexible part in the same order as the variables */
4652  int* flexenergies, /**< array of flexible energies in the same order as the variables */
4653  int* perm, /**< permutation of the variables w.r.t. the non-decreasing order of the earliest start times */
4654  int* ests, /**< array with earliest strart times sorted in non-decreasing order */
4655  int* lcts, /**< array with latest completion times sorted in non-decreasing order */
4656  int* coreEnergyAfterEst, /**< core energy after the earliest start times */
4657  int* coreEnergyAfterLct, /**< core energy after the latest completion times */
4658  SCIP_Bool* initialized, /**< was conflict analysis initialized */
4659  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
4660  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
4661  )
4662 {
4663  int coreEnergyAfterEnd;
4664  SCIP_Longint maxavailable;
4665  SCIP_Longint minavailable;
4666  SCIP_Longint totalenergy;
4667  int nests;
4668  int est;
4669  int lct;
4670  int start;
4671  int end;
4672  int v;
4673
4674  est = INT_MAX;
4675  lct = INT_MIN;
4676
4677  /* compute earliest start and latest completion time of all jobs */
4678  for( v = 0; v < nvars; ++v )
4679  {
4680  start = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[v]));
4681  end = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[v])) + durations[v];
4682
4683  est = MIN(est, start);
4684  lct = MAX(lct, end);
4685  }
4686
4687  /* adjust the effective time horizon */
4688  hmin = MAX(hmin, est);
4689  hmax = MIN(hmax, lct);
4690
4691  end = hmax + 1;
4692  coreEnergyAfterEnd = -1;
4693
4694  maxavailable = ((SCIP_Longint) hmax - hmin) * capacity;
4695  minavailable = maxavailable;
4696  totalenergy = computeTotalEnergy(durations, demands, nvars);
4697
4698  /* check if the smallest interval has a size such that the total energy fits, if so we can skip the propagator */
4699  if( ((SCIP_Longint) lcts[0] - ests[nvars-1]) * capacity >= totalenergy )
4700  return SCIP_OKAY;
4701
4702  nests = nvars;
4703
4704  /* loop over all variable in non-increasing order w.r.t. the latest completion time; thereby, the latest completion
4705  * times define the end of the time interval under investigation
4706  */
4707  for( v = nvars-1; v >= 0 && !(*cutoff); --v )
4708  {
4709  int flexenergy;
4710  int minbegin;
4711  int lbenergy;
4712  int lbcand;
4713  int i;
4714
4715  lct = lcts[v];
4716
4717  /* if the latest completion time is larger then hmax an infeasibility cannot be detected, since after hmax an
4718  * infinity capacity is available; hence we skip that
4719  */
4720  if( lct > hmax )
4721  continue;
4722
4723  /* if the latest completion time is smaller then hmin we have to stop */
4724  if( lct <= hmin )
4725  {
4726  assert(v == 0 || lcts[v-1] <= lcts[v]);
4727  break;
4728  }
4729
4730  /* if the latest completion time equals to previous end time, we can continue since this particular interval
4731  * induced by end was just analyzed
4732  */
4733  if( lct == end )
4734  continue;
4735
4736  assert(lct < end);
4737
4738  /* In case we only want to detect an overload (meaning no bound propagation) we can skip the interval; this is
4739  * the case if the free energy (the energy which is not occupied by any core) is smaller than the previous minimum
4740  * free energy; if so it means that in the next iterate the free-energy cannot be negative
4741  */
4742  if( !conshdlrdata->ttefinfer && end <= hmax && minavailable < maxavailable )
4743  {
4744  SCIP_Longint freeenergy;
4745
4746  assert(coreEnergyAfterLct[v] >= coreEnergyAfterEnd);
4747  assert(coreEnergyAfterEnd >= 0);
4748
4749  /* compute the energy which is not consumed by the cores with in the interval [lct, end) */
4750  freeenergy = capacity * ((SCIP_Longint) end - lct) - coreEnergyAfterLct[v] + coreEnergyAfterEnd;
4751
4752  if( freeenergy <= minavailable )
4753  {
4754  SCIPdebugMsg(scip, "skip latest completion time <%d> (minimum available energy <%d>, free energy <%d>)\n", lct, minavailable, freeenergy);
4755  continue;
4756  }
4757  }
4758
4759  SCIPdebugMsg(scip, "check intervals ending with <%d>\n", lct);
4760
4761  end = lct;
4762  coreEnergyAfterEnd = coreEnergyAfterLct[v];
4763
4764  flexenergy = 0;
4765  minavailable = maxavailable;
4766  minbegin = hmax;
4767  lbcand = -1;
4768  lbenergy = 0;
4769
4770  /* loop over the job in non-increasing order w.r.t. the earliest start time; these earliest start time are
4771  * defining the beginning of the time interval under investigation; Thereby, the time interval gets wider and
4772  * wider
4773  */
4774  for( i = nests-1; i >= 0; --i )
4775  {
4776  SCIP_VAR* var;
4777  SCIP_Longint freeenergy;
4778  int duration;
4779  int demand;
4780  int begin;
4781  int idx;
4782  int lst;
4783
4784  idx = perm[i];
4785  assert(idx >= 0);
4786  assert(idx < nvars);
4787  assert(!(*cutoff));
4788
4789  /* the earliest start time of the job */
4790  est = ests[i];
4791
4792  /* if the job starts after the current end, we can skip it and do not need to consider it again since the
4793  * latest completion times (which define end) are scant in non-increasing order
4794  */
4795  if( end <= est )
4796  {
4797  nests--;
4798  continue;
4799  }
4800
4801  /* check if the interval has a size such that the total energy fits, if so we can skip all intervals with the
4802  * current ending time
4803  */
4804  if( ((SCIP_Longint) end - est) * capacity >= totalenergy )
4805  break;
4806
4807  var = vars[idx];
4808  assert(var != NULL);
4809
4810  duration = durations[idx];
4811  assert(duration > 0);
4812
4813  demand = demands[idx];
4814  assert(demand > 0);
4815
4816  lct = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var)) + duration;
4817
4818  /* the latest start time of the free part of the job */
4819  lst = lsts[idx];
4820
4821  /* in case the earliest start time is equal to minbegin, the job lies completely within the time window under
4822  * investigation; hence the overload check will do the the job
4823  */
4824  assert(est <= minbegin);
4825  if( minavailable < maxavailable && est < minbegin )
4826  {
4827  assert(!(*cutoff));
4828
4829  /* try to tighten the upper bound */
4830  SCIP_CALL( tightenUbTTEF(scip, conshdlrdata, nvars, vars, durations, demands, capacity, hmin, hmax,
4831  var, duration, demand, est, lst, lct, minbegin, end, minavailable, &(newubs[idx]), &(ubinferinfos[idx]),
4832  initialized, explanation, cutoff) );
4833
4834  if( *cutoff )
4835  break;
4836  }
4837
4838  SCIPdebugMsg(scip, "check variable <%s>[%g,%g] (duration %d, demands %d, est <%d>, lst of free part <%d>\n",
4839  SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), duration, demand, est, lst);
4840
4841  begin = est;
4842  assert(SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var)) == est);
4843
4844  /* if the earliest start time is smaller than hmin we can stop here since the next job will not decrease the
4845  * free energy
4846  */
4847  if( begin < hmin )
4848  break;
4849
4850  /* compute the contribution to the flexible energy */
4851  if( lct <= end )
4852  {
4853  /* if the jobs has to finish before the end, all the energy has to be scheduled */
4854  assert(lst >= begin);
4855  assert(flexenergies[idx] >= 0);
4856  flexenergy += flexenergies[idx];
4857  }
4858  else
4859  {
4860  /* the job partly overlaps with the end */
4861  int candenergy;
4862  int energy;
4863
4864  /* compute the flexible energy which is part of the time interval for sure if the job is scheduled
4865  * w.r.t. latest start time
4866  *
4867  * @note we need to be aware of the effective horizon
4868  */
4869  energy = MIN(flexenergies[idx], demands[idx] * MAX(0, (end - lst)));
4870  assert(end - lst < duration);
4871  assert(energy >= 0);
4872
4873  /* adjust the flexible energy of the time interval */
4874  flexenergy += energy;
4875
4876  /* compute the flexible energy of the job which is not part of flexible energy of the time interval */
4877  candenergy = MIN(flexenergies[idx], demands[idx] * (end - begin)) - energy;
4878  assert(candenergy >= 0);
4879
4880  /* check if we found a better candidate */
4881  if( candenergy > lbenergy )
4882  {
4883  lbenergy = candenergy;
4884  lbcand = idx;
4885  }
4886  }
4887
4888  SCIPdebugMsg(scip, "time window [%d,%d) flexible energy <%d>\n", begin, end, flexenergy);
4889  assert(coreEnergyAfterEst[i] >= coreEnergyAfterEnd);
4890
4891  /* compute the energy which is not used yet */
4892  freeenergy = capacity * ((SCIP_Longint) end - begin) - flexenergy - coreEnergyAfterEst[i] + coreEnergyAfterEnd;
4893
4894  /* check overload */
4895  if( freeenergy < 0 )
4896  {
4897  SCIPdebugMsg(scip, "analyze overload within time window [%d,%d) capacity %d\n", begin, end, capacity);
4898
4899  /* initialize conflict analysis if conflict analysis is applicable */
4901  {
4902  /* analyze infeasibilty */
4904
4905  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
4906  begin, end, NULL, SCIP_BOUNDTYPE_UPPER, NULL, SCIP_UNKNOWN,
4907  conshdlrdata->usebdwidening, explanation) );
4908
4909  (*initialized) = TRUE;
4910  }
4911
4912  (*cutoff) = TRUE;
4913
4914  /* for the statistic we count the number of times a cutoff was detected due the time-time-edge-finding */
4915  SCIPstatistic( SCIPconshdlrGetData(SCIPfindConshdlr(scip, CONSHDLR_NAME))->ncutoffoverloadTTEF++ );
4916
4917  break;
4918  }
4919
4920  /* check if the available energy is not sufficent to schedule the flexible energy of the best candidate job */
4921  if( lbenergy > 0 && freeenergy < lbenergy )
4922  {
4923  SCIP_Longint energy;
4924  int newlb;
4925  int ect;
4926
4927  ect = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[lbcand])) + durations[lbcand];
4928  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[lbcand]));
4929
4930  /* remove the energy of our job from the ... */
4931  energy = freeenergy + (computeCoreWithInterval(begin, end, ect, lst) + MAX(0, (SCIP_Longint) end - lsts[lbcand])) * demands[lbcand];
4932
4933  newlb = end - (int)(energy / demands[lbcand]);
4934
4935  if( newlb > lst )
4936  {
4937  /* initialize conflict analysis if conflict analysis is applicable */
4939  {
4940  SCIP_Real relaxedbd;
4941
4942  /* analyze infeasibilty */
4944
4945  relaxedbd = lst + 1.0;
4946
4947  /* added to upper bound (which was overcut be new lower bound) of the variable */
4948  SCIP_CALL( SCIPaddConflictUb(scip, vars[lbcand], NULL) );
4949
4950  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
4951  begin, end, vars[lbcand], SCIP_BOUNDTYPE_LOWER, NULL, relaxedbd,
4952  conshdlrdata->usebdwidening, explanation) );
4953
4954  (*initialized) = TRUE;
4955  }
4956
4957  (*cutoff) = TRUE;
4958  break;
4959  }
4960  else if( newlb > newlbs[lbcand] )
4961  {
4962  INFERINFO inferinfo;
4963
4964  /* construct inference information */
4965  inferinfo = getInferInfo(PROPRULE_3_TTEF, begin, end);
4966
4967  /* buffer upper bound change */
4968  lbinferinfos[lbcand] = inferInfoToInt(inferinfo);
4969  newlbs[lbcand] = newlb;
4970  }
4971  }
4972
4973  /* check if the current interval has a smaller free energy */
4974  if( minavailable > freeenergy )
4975  {
4976  minavailable = freeenergy;
4977  minbegin = begin;
4978  }
4979  assert(minavailable >= 0);
4980  }
4981  }
4982
4983  return SCIP_OKAY;
4984 }
4985
4986 /** propagate the lower bounds and "opportunistically" the upper bounds using the time-table edge-finding algorithm */
4987 static
4989  SCIP* scip, /**< SCIP data structure */
4990  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
4991  int nvars, /**< number of start time variables (activities) */
4992  SCIP_VAR** vars, /**< array of start time variables */
4993  int* durations, /**< array of durations */
4994  int* demands, /**< array of demands */
4995  int capacity, /**< cumulative capacity */
4996  int hmin, /**< left bound of time axis to be considered (including hmin) */
4997  int hmax, /**< right bound of time axis to be considered (not including hmax) */
4998  int* newlbs, /**< array to buffer new lower bounds */
4999  int* newubs, /**< array to buffer new upper bounds */
5000  int* lbinferinfos, /**< array to store the inference information for the lower bound changes */
5001  int* ubinferinfos, /**< array to store the inference information for the upper bound changes */
5002  int* ects, /**< array of earliest completion time of the flexible part in the same order as the variables */
5003  int* flexenergies, /**< array of flexible energies in the same order as the variables */
5004  int* perm, /**< permutation of the variables w.r.t. the non-decreasing order of the latest completion times */
5005  int* ests, /**< array with earliest strart times sorted in non-decreasing order */
5006  int* lcts, /**< array with latest completion times sorted in non-decreasing order */
5007  int* coreEnergyAfterEst, /**< core energy after the earliest start times */
5008  int* coreEnergyAfterLct, /**< core energy after the latest completion times */
5009  SCIP_Bool* initialized, /**< was conflict analysis initialized */
5010  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
5011  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
5012  )
5013 {
5014  int coreEnergyAfterStart;
5015  SCIP_Longint maxavailable;
5016  SCIP_Longint minavailable;
5017  SCIP_Longint totalenergy;
5018  int nlcts;
5019  int begin;
5020  int minest;
5021  int maxlct;
5022  int start;
5023  int end;
5024  int v;
5025
5026  if( *cutoff )
5027  return SCIP_OKAY;
5028
5029  begin = hmin - 1;
5030
5031  minest = INT_MAX;
5032  maxlct = INT_MIN;
5033
5034  /* compute earliest start and latest completion time of all jobs */
5035  for( v = 0; v < nvars; ++v )
5036  {
5037  start = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[v]));
5038  end = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[v])) + durations[v];
5039
5040  minest = MIN(minest, start);
5041  maxlct = MAX(maxlct, end);
5042  }
5043
5044  /* adjust the effective time horizon */
5045  hmin = MAX(hmin, minest);
5046  hmax = MIN(hmax, maxlct);
5047
5048  maxavailable = ((SCIP_Longint) hmax - hmin) * capacity;
5049  totalenergy = computeTotalEnergy(durations, demands, nvars);
5050
5051  /* check if the smallest interval has a size such that the total energy fits, if so we can skip the propagator */
5052  if( ((SCIP_Longint) lcts[0] - ests[nvars-1]) * capacity >= totalenergy )
5053  return SCIP_OKAY;
5054
5055  nlcts = 0;
5056
5057  /* loop over all variable in non-decreasing order w.r.t. the earliest start times; thereby, the earliest start times
5058  * define the start of the time interval under investigation
5059  */
5060  for( v = 0; v < nvars; ++v )
5061  {
5062  int flexenergy;
5063  int minend;
5064  int ubenergy;
5065  int ubcand;
5066  int est;
5067  int i;
5068
5069  est = ests[v];
5070
5071  /* if the earliest start time is smaller then hmin an infeasibility cannot be detected, since before hmin an
5072  * infinity capacity is available; hence we skip that
5073  */
5074  if( est < hmin )
5075  continue;
5076
5077  /* if the earliest start time is larger or equal then hmax we have to stop */
5078  if( est >= hmax )
5079  break;
5080
5081  /* if the latest earliest start time equals to previous start time, we can continue since this particular interval
5082  * induced by start was just analyzed
5083  */
5084  if( est == begin )
5085  continue;
5086
5087  assert(est > begin);
5088
5089  SCIPdebugMsg(scip, "check intervals starting with <%d>\n", est);
5090
5091  begin = est;
5092  coreEnergyAfterStart = coreEnergyAfterEst[v];
5093
5094  flexenergy = 0;
5095  minavailable = maxavailable;
5096  minend = hmin;
5097  ubcand = -1;
5098  ubenergy = 0;
5099
5100  /* loop over the job in non-decreasing order w.r.t. the latest completion time; these latest completion times are
5101  * defining the ending of the time interval under investigation; thereby, the time interval gets wider and wider
5102  */
5103  for( i = nlcts; i < nvars; ++i )
5104  {
5105  SCIP_VAR* var;
5106  SCIP_Longint freeenergy;
5107  int duration;
5108  int demand;
5109  int idx;
5110  int lct;
5111  int ect;
5112
5113  idx = perm[i];
5114  assert(idx >= 0);
5115  assert(idx < nvars);
5116  assert(!(*cutoff));
5117
5118  /* the earliest start time of the job */
5119  lct = lcts[i];
5120
5121  /* if the job has a latest completion time before the the current start, we can skip it and do not need to
5122  * consider it again since the earliest start times (which define the start) are scant in non-decreasing order
5123  */
5124  if( lct <= begin )
5125  {
5126  nlcts++;
5127  continue;
5128  }
5129
5130  /* check if the interval has a size such that the total energy fits, if so we can skip all intervals which
5131  * start with current beginning time
5132  */
5133  if( ((SCIP_Longint) lct - begin) * capacity >= totalenergy )
5134  break;
5135
5136  var = vars[idx];
5137  assert(var != NULL);
5138
5139  duration = durations[idx];
5140  assert(duration > 0);
5141
5142  demand = demands[idx];
5143  assert(demand > 0);
5144
5145  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
5146
5147  /* the earliest completion time of the flexible part of the job */
5148  ect = ects[idx];
5149
5150  /* in case the latest completion time is equal to minend, the job lies completely within the time window under
5151  * investigation; hence the overload check will do the the job
5152  */
5153  assert(lct >= minend);
5154  if( minavailable < maxavailable && lct > minend )
5155  {
5156  assert(!(*cutoff));
5157
5158  /* try to tighten the upper bound */
5159  SCIP_CALL( tightenLbTTEF(scip, conshdlrdata, nvars, vars, durations, demands, capacity, hmin, hmax,
5160  var, duration, demand, est, ect, lct, begin, minend, minavailable, &(newlbs[idx]), &(lbinferinfos[idx]),
5161  initialized, explanation, cutoff) );
5162
5163  if( *cutoff )
5164  return SCIP_OKAY;
5165  }
5166
5167  SCIPdebugMsg(scip, "check variable <%s>[%g,%g] (duration %d, demands %d, est <%d>, ect of free part <%d>\n",
5168  SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), duration, demand, est, ect);
5169
5170  end = lct;
5171  assert(SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var)) + duration == lct);
5172
5173  /* if the latest completion time is larger than hmax we can stop here since the next job will not decrease the
5174  * free energy
5175  */
5176  if( end > hmax )
5177  break;
5178
5179  /* compute the contribution to the flexible energy */
5180  if( est >= begin )
5181  {
5182  /* if the jobs has to finish before the end, all the energy has to be scheduled */
5183  assert(ect <= end);
5184  assert(flexenergies[idx] >= 0);
5185  flexenergy += flexenergies[idx];
5186  }
5187  else
5188  {
5189  /* the job partly overlaps with the end */
5190  int candenergy;
5191  int energy;
5192
5193  /* compute the flexible energy which is part of the time interval for sure if the job is scheduled
5194  * w.r.t. latest start time
5195  *
5196  * @note we need to be aware of the effective horizon
5197  */
5198  energy = MIN(flexenergies[idx], demands[idx] * MAX(0, (ect - begin)));
5199  assert(ect - begin < duration);
5200  assert(energy >= 0);
5201
5202  /* adjust the flexible energy of the time interval */
5203  flexenergy += energy;
5204
5205  /* compute the flexible energy of the job which is not part of flexible energy of the time interval */
5206  candenergy = MIN(flexenergies[idx], demands[idx] * (end - begin)) - energy;
5207  assert(candenergy >= 0);
5208
5209  /* check if we found a better candidate */
5210  if( candenergy > ubenergy )
5211  {
5212  ubenergy = candenergy;
5213  ubcand = idx;
5214  }
5215  }
5216
5217  SCIPdebugMsg(scip, "time window [%d,%d) flexible energy <%d>\n", begin, end, flexenergy);
5218  assert(coreEnergyAfterLct[i] <= coreEnergyAfterStart);
5219
5220  /* compute the energy which is not used yet */
5221  freeenergy = capacity * ((SCIP_Longint) end - begin) - flexenergy - coreEnergyAfterStart + coreEnergyAfterLct[i];
5222
5223  /* check overload */
5224  if( freeenergy < 0 )
5225  {
5226  SCIPdebugMsg(scip, "analyze overload within time window [%d,%d) capacity %d\n", begin, end, capacity);
5227
5228  /* initialize conflict analysis if conflict analysis is applicable */
5230  {
5231  /* analyze infeasibilty */
5233
5234  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
5235  begin, end, NULL, SCIP_BOUNDTYPE_UPPER, NULL, SCIP_UNKNOWN,
5236  conshdlrdata->usebdwidening, explanation) );
5237
5238  (*initialized) = TRUE;
5239  }
5240
5241  (*cutoff) = TRUE;
5242
5243  /* for the statistic we count the number of times a cutoff was detected due the time-time-edge-finding */
5244  SCIPstatistic( SCIPconshdlrGetData(SCIPfindConshdlr(scip, CONSHDLR_NAME))->ncutoffoverloadTTEF++ );
5245
5246  return SCIP_OKAY;
5247  }
5248
5249  /* check if the available energy is not sufficent to schedule the flexible energy of the best candidate job */
5250  if( ubenergy > 0 && freeenergy < ubenergy )
5251  {
5252  SCIP_Longint energy;
5253  int newub;
5254  int lst;
5255
5256  duration = durations[ubcand];
5257  assert(duration > 0);
5258
5259  ect = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[ubcand])) + duration;
5260  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[ubcand]));
5261
5262  /* remove the energy of our job from the ... */
5263  energy = freeenergy + (computeCoreWithInterval(begin, end, ect, lst) + MAX(0, (SCIP_Longint) ects[ubcand] - begin)) * demands[ubcand];
5264
5265  newub = begin - duration + (int)(energy / demands[ubcand]);
5266
5267  if( newub < ect - duration )
5268  {
5269  /* initialize conflict analysis if conflict analysis is applicable */
5271  {
5272  SCIP_Real relaxedbd;
5273  /* analyze infeasibilty */
5275
5276  relaxedbd = ect - duration - 1.0;
5277
5278  /* added to lower bound (which was undercut be new upper bound) of the variable */
5279  SCIP_CALL( SCIPaddConflictUb(scip, vars[ubcand], NULL) );
5280
5281  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
5282  begin, end, vars[ubcand], SCIP_BOUNDTYPE_UPPER, NULL, relaxedbd,
5283  conshdlrdata->usebdwidening, explanation) );
5284
5285  (*initialized) = TRUE;
5286  }
5287
5288  (*cutoff) = TRUE;
5289  return SCIP_OKAY;
5290  }
5291  else if( newub < newubs[ubcand] )
5292  {
5293  INFERINFO inferinfo;
5294
5295  /* construct inference information */
5296  inferinfo = getInferInfo(PROPRULE_3_TTEF, begin, end);
5297
5298  /* buffer upper bound change */
5299  ubinferinfos[ubcand] = inferInfoToInt(inferinfo);
5300  newubs[ubcand] = newub;
5301  }
5302  }
5303
5304  /* check if the current interval has a smaller free energy */
5305  if( minavailable > freeenergy )
5306  {
5307  minavailable = freeenergy;
5308  minend = end;
5309  }
5310  assert(minavailable >= 0);
5311  }
5312  }
5313
5314  return SCIP_OKAY;
5315 }
5316
5317 /** checks whether the instance is infeasible due to a overload within a certain time frame using the idea of time-table
5318  * edge-finding
5319  *
5320  * @note The algorithm is based on the following two papers:
5321  * - Petr Vilim, "Timetable Edge Finding Filtering Algorithm for Discrete Cumulative Resources", In: Tobias
5322  * Achterberg and J. Christopher Beck (Eds.), Integration of AI and OR Techniques in Constraint Programming for
5323  * Combinatorial Optimization Problems (CPAIOR 2011), LNCS 6697, pp 230--245
5324  * - Andreas Schutt, Thibaut Feydy, and Peter J. Stuckey, "Explaining Time-Table-Edge-Finding Propagation for the
5325  * Cumulative Resource Constraint (submitted to CPAIOR 2013)
5326  */
5327 static
5329  SCIP* scip, /**< SCIP data structure */
5330  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
5331  SCIP_PROFILE* profile, /**< current core profile */
5332  int nvars, /**< number of start time variables (activities) */
5333  SCIP_VAR** vars, /**< array of start time variables */
5334  int* durations, /**< array of durations */
5335  int* demands, /**< array of demands */
5336  int capacity, /**< cumulative capacity */
5337  int hmin, /**< left bound of time axis to be considered (including hmin) */
5338  int hmax, /**< right bound of time axis to be considered (not including hmax) */
5339  SCIP_CONS* cons, /**< constraint which is propagated (needed to SCIPinferVar**Cons()) */
5340  int* nchgbds, /**< pointer to store the number of bound changes */
5341  SCIP_Bool* initialized, /**< was conflict analysis initialized */
5342  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
5343  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
5344  )
5345 {
5346  int* coreEnergyAfterEst;
5347  int* coreEnergyAfterLct;
5348  int* flexenergies;
5349  int* permests;
5350  int* permlcts;
5351  int* lcts;
5352  int* ests;
5353  int* ects;
5354  int* lsts;
5355
5356  int* newlbs;
5357  int* newubs;
5358  int* lbinferinfos;
5359  int* ubinferinfos;
5360
5361  int v;
5362
5363  /* check if a cutoff was already detected */
5364  if( (*cutoff) )
5365  return SCIP_OKAY;
5366
5367  /* check if at least the basic overload checking should be perfomed */
5368  if( !conshdlrdata->ttefcheck )
5369  return SCIP_OKAY;
5370
5371  SCIPdebugMsg(scip, "run time-table edge-finding overload checking\n");
5372
5373  SCIP_CALL( SCIPallocBufferArray(scip, &coreEnergyAfterEst, nvars) );
5374  SCIP_CALL( SCIPallocBufferArray(scip, &coreEnergyAfterLct, nvars) );
5375  SCIP_CALL( SCIPallocBufferArray(scip, &flexenergies, nvars) );
5376  SCIP_CALL( SCIPallocBufferArray(scip, &permlcts, nvars) );
5377  SCIP_CALL( SCIPallocBufferArray(scip, &permests, nvars) );
5378  SCIP_CALL( SCIPallocBufferArray(scip, &lcts, nvars) );
5379  SCIP_CALL( SCIPallocBufferArray(scip, &ests, nvars) );
5380  SCIP_CALL( SCIPallocBufferArray(scip, &ects, nvars) );
5381  SCIP_CALL( SCIPallocBufferArray(scip, &lsts, nvars) );
5382
5383  SCIP_CALL( SCIPallocBufferArray(scip, &newlbs, nvars) );
5384  SCIP_CALL( SCIPallocBufferArray(scip, &newubs, nvars) );
5385  SCIP_CALL( SCIPallocBufferArray(scip, &lbinferinfos, nvars) );
5386  SCIP_CALL( SCIPallocBufferArray(scip, &ubinferinfos, nvars) );
5387
5388  /* we need to buffer the bound changes since the propagation algorithm cannot handle new bound dynamically */
5389  for( v = 0; v < nvars; ++v )
5390  {
5391  newlbs[v] = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(vars[v]));
5392  newubs[v] = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(vars[v]));
5393  lbinferinfos[v] = 0;
5394  ubinferinfos[v] = 0;
5395  }
5396
5397  /* collect earliest start times, latest completion time, and free energy contributions */
5398  collectDataTTEF(scip, nvars, vars, durations, demands, hmin, hmax, permests, ests, permlcts, lcts, ects, lsts, flexenergies);
5399
5400  /* sort the earliest start times and latest completion in non-decreasing order */
5401  SCIPsortIntInt(ests, permests, nvars);
5402  SCIPsortIntInt(lcts, permlcts, nvars);
5403
5404  /* compute for the different earliest start and latest completion time the core energy of the corresponding time
5405  * points
5406  */
5407  SCIP_CALL( computeCoreEngeryAfter(scip, profile, nvars, ests, lcts, coreEnergyAfterEst, coreEnergyAfterLct) );
5408
5409  /* propagate the upper bounds and "opportunistically" the lower bounds */
5410  SCIP_CALL( propagateUbTTEF(scip, conshdlrdata, nvars, vars, durations, demands, capacity, hmin, hmax,
5411  newlbs, newubs, lbinferinfos, ubinferinfos, lsts, flexenergies,
5412  permests, ests, lcts, coreEnergyAfterEst, coreEnergyAfterLct, initialized, explanation, cutoff) );
5413
5414  /* propagate the lower bounds and "opportunistically" the upper bounds */
5415  SCIP_CALL( propagateLbTTEF(scip, conshdlrdata, nvars, vars, durations, demands, capacity, hmin, hmax,
5416  newlbs, newubs, lbinferinfos, ubinferinfos, ects, flexenergies,
5417  permlcts, ests, lcts, coreEnergyAfterEst, coreEnergyAfterLct, initialized, explanation, cutoff) );
5418
5419  /* apply the buffer bound changes */
5420  for( v = 0; v < nvars && !(*cutoff); ++v )
5421  {
5422  SCIP_Bool infeasible;
5423  SCIP_Bool tightened;
5424
5425  SCIP_CALL( SCIPinferVarLbCons(scip, vars[v], (SCIP_Real)newlbs[v], cons, lbinferinfos[v], TRUE, &infeasible, &tightened) );
5426
5427  /* since we change first the lower bound of the variable an infeasibilty should be detected */
5428  assert(!infeasible);
5429
5430  if( tightened )
5431  {
5432  (*nchgbds)++;
5433
5434  /* for the statistic we count the number of times a cutoff was detected due the time-time */
5436  }
5437
5438  SCIP_CALL( SCIPinferVarUbCons(scip, vars[v], (SCIP_Real)newubs[v], cons, ubinferinfos[v], TRUE, &infeasible, &tightened) );
5439
5440  /* since upper bound was compute w.r.t. the "old" bound the previous lower bound update together with this upper
5441  * bound update can be infeasible
5442  */
5443  if( infeasible )
5444  {
5446  {
5447  INFERINFO inferinfo;
5448  SCIP_VAR* var;
5449  int begin;
5450  int end;
5451
5452  var = vars[v];
5453  assert(var != NULL);
5454
5455  /* initialize conflict analysis */
5457
5458  /* convert int to inference information */
5459  inferinfo = intToInferInfo(ubinferinfos[v]);
5460
5461  /* collect time window from inference information */
5462  begin = inferInfoGetData1(inferinfo);
5463  end = inferInfoGetData2(inferinfo);
5464  assert(begin < end);
5465
5466  /* added to lower bound (which was undercut be new upper bound) of the variable */
5467  SCIP_CALL( SCIPaddConflictLb(scip, var, NULL) );
5468
5469  /* analysis the upper bound change */
5470  SCIP_CALL( analyzeEnergyRequirement(scip, nvars, vars, durations, demands, capacity,
5471  begin, end, var, SCIP_BOUNDTYPE_UPPER, NULL, SCIPvarGetLbLocal(vars[v]) - 1.0,
5472  conshdlrdata->usebdwidening, explanation) );
5473
5474  (*initialized) = TRUE;
5475  }
5476
5477  /* for the statistic we count the number of times a cutoff was detected due the time-time */
5478  SCIPstatistic( SCIPconshdlrGetData(SCIPfindConshdlr(scip, CONSHDLR_NAME))->ncutoffoverloadTTEF++ );
5479
5480  (*cutoff) = TRUE;
5481  break;
5482  }
5483
5484  if( tightened )
5485  {
5486  (*nchgbds)++;
5487
5488  /* for the statistic we count the number of times a cutoff was detected due the time-time */
5490  }
5491  }
5492
5493  SCIPfreeBufferArray(scip, &ubinferinfos);
5494  SCIPfreeBufferArray(scip, &lbinferinfos);
5495  SCIPfreeBufferArray(scip, &newubs);
5496  SCIPfreeBufferArray(scip, &newlbs);
5497
5498  /* free buffer arrays */
5499  SCIPfreeBufferArray(scip, &lsts);
5500  SCIPfreeBufferArray(scip, &ects);
5501  SCIPfreeBufferArray(scip, &ests);
5502  SCIPfreeBufferArray(scip, &lcts);
5503  SCIPfreeBufferArray(scip, &permests);
5504  SCIPfreeBufferArray(scip, &permlcts);
5505  SCIPfreeBufferArray(scip, &flexenergies);
5506  SCIPfreeBufferArray(scip, &coreEnergyAfterLct);
5507  SCIPfreeBufferArray(scip, &coreEnergyAfterEst);
5508
5509  return SCIP_OKAY;
5510 }
5511
5512 /** a cumulative condition is not satisfied if its capacity is exceeded at a time where jobs cannot be shifted (core)
5513  * anymore we build up a cumulative profile of all cores of jobs and try to improve bounds of all jobs; also known as
5514  * time table propagator
5515  */
5516 static
5518  SCIP* scip, /**< SCIP data structure */
5519  SCIP_CONSHDLRDATA* conshdlrdata, /**< constraint handler data */
5520  SCIP_PROFILE* profile, /**< core profile */
5521  int nvars, /**< number of start time variables (activities) */
5522  SCIP_VAR** vars, /**< array of start time variables */
5523  int* durations, /**< array of durations */
5524  int* demands, /**< array of demands */
5525  int capacity, /**< cumulative capacity */
5526  int hmin, /**< left bound of time axis to be considered (including hmin) */
5527  int hmax, /**< right bound of time axis to be considered (not including hmax) */
5528  SCIP_CONS* cons, /**< constraint which is propagated (needed to SCIPinferVar**Cons()) */
5529  int* nchgbds, /**< pointer to store the number of bound changes */
5530  SCIP_Bool* initialized, /**< was conflict analysis initialized */
5531  SCIP_Bool* explanation, /**< bool array which marks the variable which are part of the explanation if a cutoff was detected, or NULL */
5532  SCIP_Bool* cutoff /**< pointer to store if the constraint is infeasible */
5533  )
5534 {
5535  SCIP_Bool infeasible;
5536  int v;
5537
5538  assert(scip != NULL);
5539  assert(nvars > 0);
5540  assert(cons != NULL);
5541  assert(cutoff != NULL);
5542
5543  /* check if already a cutoff was detected */
5544  if( (*cutoff) )
5545  return SCIP_OKAY;
5546
5547  /* check if the time tabling should infer bounds */
5548  if( !conshdlrdata->ttinfer )
5549  return SCIP_OKAY;
5550
5551  assert(*initialized == FALSE);
5552
5553  SCIPdebugMsg(scip, "propagate cores of cumulative condition of constraint <%s>[%d,%d) <= %d\n",
5554  SCIPconsGetName(cons), hmin, hmax, capacity);
5555
5556  infeasible = FALSE;
5557
5558  /* if core profile is empty; nothing to do */
5559  if( SCIPprofileGetNTimepoints(profile) <= 1 )
5560  return SCIP_OKAY;
5561
5562  /* start checking each job whether the bounds can be improved */
5563  for( v = 0; v < nvars; ++v )
5564  {
5565  SCIP_VAR* var;
5566  int demand;
5567  int duration;
5568  int begin;
5569  int end;
5570  int est;
5571  int lst;
5572
5573  var = vars[v];
5574  assert(var != NULL);
5575
5576  duration = durations[v];
5577  assert(duration > 0);
5578
5579  /* collect earliest and latest start time */
5580  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
5581  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
5582
5583  /* check if the start time variables is already fixed; in that case we can ignore the job */
5584  if( est == lst )
5585  continue;
5586
5587  /* check if the job runs completely outside of the effective horizon [hmin, hmax); if so skip it */
5588  if( lst + duration <= hmin || est >= hmax )
5589  continue;
5590
5591  /* compute core interval w.r.t. effective time horizon */
5592  begin = MAX(hmin, lst);
5593  end = MIN(hmax, est + duration);
5594
5595  demand = demands[v];
5596  assert(demand > 0);
5597
5598  /* if the job has a core, remove it first */
5599  if( begin < end )
5600  {
5601  SCIPdebugMsg(scip, "variable <%s>[%g,%g] (duration %d, demand %d): remove core [%d,%d)\n",
5602  SCIPvarGetName(var), SCIPvarGetLbLocal(var), SCIPvarGetUbLocal(var), duration, demand, begin, end);
5603
5604  SCIP_CALL( SCIPprofileDeleteCore(profile, begin, end, demand) );
5605  }
5606
5607  /* first try to update the earliest start time */
5608  SCIP_CALL( coretimesUpdateLb(scip, nvars, vars, durations, demands, capacity, hmin, hmax, cons,
5609  profile, v, nchgbds, conshdlrdata->usebdwidening, initialized, explanation, cutoff) );
5610
5611  if( *cutoff )
5612  break;
5613
5614  /* second try to update the latest start time */
5615  SCIP_CALL( coretimesUpdateUb(scip, var, duration, demand, capacity, cons,
5616  profile, v, nchgbds) );
5617
5618  if( *cutoff )
5619  break;
5620
5621  /* collect the potentially updated earliest and latest start time */
5622  est = SCIPconvertRealToInt(scip, SCIPvarGetLbLocal(var));
5623  lst = SCIPconvertRealToInt(scip, SCIPvarGetUbLocal(var));
5624
5625  /* compute core interval w.r.t. effective time horizon */
5626  begin = MAX(hmin, lst);
5627  end = MIN(hmax, est + duration);
5628
5629  /* after updating the bound we might have a new core */
5630  if( begin < end )
5631  {
5632  int pos;
5633
5634  SCIPdebugMsg(scip, "variable <%s>[%d,%d] (duration %d, demand %d): add core [%d,%d)\n",
5635  SCIPvarGetName(var), est, lst, duration, demand, begin, end);
5636
5637  SCIP_CALL( SCIPprofileInsertCore(profile, begin, end, demand, &pos, &infeasible) );
5638
5639  if( infeasible )
5640  {
5641  /* use conflict analysis to analysis the core insertion which was infeasible */
5642  SCIP_CALL( analyseInfeasibelCoreInsertion(scip, nvars, vars, durations, demands, capacity, hmin, hmax,
5643  var, duration, demand, SCIPprofileGetTime(profile, pos), conshdlrdata->usebdwidening, initialized, explanation) );
5644
5645  if( explanation != NULL )
5646  explanation[v] = TRUE;
5647
5648  (*cutoff) = TRUE;
5649
5650  /* for the statistic we count the number of times a cutoff was detected due the time-time */
5651  SCIPstatistic( SCIPconshdlrGetData(SCIPfindConshdlr(scip, CONSHDLR_NAME))->ncutofftimetable++ );
5652
5653  break;
5654  }
5655  }
5656  }
5657
5658  return SCIP_OKAY;
5659 }
5660
5661
5662 /** node data structure for the binary tree used for edgefinding (with overload checking) */
5663 struct SCIP_NodeData
5664 {
5665  SCIP_VAR* var; /**< start time variable of the job if the node data belongs to a leaf, otherwise NULL */
5666  SCIP_Real key; /**< key which is to insert the corresponding search node */
5667  int est; /**< earliest start time if the node data belongs to a leaf */
5668  int lct; /**< latest completion time if the node data belongs to a leaf */
5669  int demand; /**< demand of the job if the node data belongs to a leaf */
5670  int duration; /**< duration of the job if the node data belongs to a leaf */
5671  int leftadjust; /**< left adjustments of the duration w.r.t. hmin */
5672  int rightadjust; /**< right adjustments of the duration w.r.t. hmax */
5673  SCIP_Longint enveloptheta; /**< the maximal energy of a subset of jobs part of the theta set */
5674  int energytheta; /**< energy of the subset of the jobs which are part of theta set */
5675  int energylambda;
5676  SCIP_Longint enveloplambda;
5677  int idx; /**< index of the start time variable in the (global) variable array */
5678  SCIP_Bool intheta; /**< belongs the node to the theta set (otherwise to the lambda set) */
5679 };
5680 typedef struct SCIP_NodeData SCIP_NODEDATA;
5682 /** creates a node data structure */
5683 static
5685  SCIP* scip, /**< SCIP data structure */
5686  SCIP_NODEDATA** nodedata /**< pointer to store the create node data */
5687  )
5688 {
5689  SCIP_CALL( SCIPallocBuffer(scip, nodedata) );
5690  (*nodedata)->var = NULL;
5691  (*nodedata)->key = SCIP_INVALID;
5692  (*nodedata)->est = INT_MIN;
5693  (*nodedata)->lct = INT_MAX;
5694  (*nodedata)->duration = 0;
5695  (*nodedata)->demand = 0;
5696  (*nodedata)->enveloptheta = -1;
5697  (*nodedata)->energytheta = 0;
5698  (*nodedata)->enveloplambda = -1;
5699  (*nodedata)->energylambda = -1;
5700  (*nodedata)->idx = -1;
5701  (*nodedata)->intheta = TRUE;
5702
5703  return SCIP_OKAY;
5704 }
5705
5706 /** frees a node data structure */
5707 static
5708 void freeNodedata(
5709  SCIP* scip, /**< SCIP data structure */
5710  SCIP_NODEDATA** nodedata /**< pointer to store node data which should be freed */
5711  )
5712 {
5713  if( *nodedata != NULL )
5714  {
5715  SCIPfreeBuffer(scip, nodedata);
5716  }
5717 }
5718
5719 /** update node data structure strating form the given node along the path to the root node */
5720 static
5722  SCIP* scip, /**< SCIP data structure */
5723  SCIP_BTNODE* node /**< search node which inserted */
5724  )
5725 {
5726  SCIP_BTNODE* left;
5727  SCIP_BTNODE* right;
5728  SCIP_NODEDATA* nodedata;
5729  SCIP_NODEDATA* leftdata;
5730  SCIP_NODEDATA* rightdata;
5731
5732  SCIPdebugMsg(scip, "update envelop starting from node <%p>\n", (void*)node);
5733
5734  if( SCIPbtnodeIsLeaf(node) )
5735  node = SCIPbtnodeGetParent(node);
5736
5737  while( node != NULL )
5738  {
5739  /* get node data */
5740  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(node);
5741  assert(nodedata != NULL);
5742
5743  /* collect node data from left node */
5744  left = SCIPbtnodeGetLeftchild(node);
5745  assert(left != NULL);
5746  leftdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(left);
5747  assert(leftdata != NULL);
5748
5749  /* collect node data from right node */
5750  right = SCIPbtnodeGetRightchild(node);
5751  assert(right != NULL);
5752  rightdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(right);
5753  assert(rightdata != NULL);
5754
5755  /* update envelop and energy */
5756  if( leftdata->enveloptheta >= 0 )
5757  {
5758  assert(rightdata->energytheta != -1);
5759  nodedata->enveloptheta = MAX(leftdata->enveloptheta + rightdata->energytheta, rightdata->enveloptheta);
5760  }
5761  else
5762  nodedata->enveloptheta = rightdata->enveloptheta;
5763
5764  assert(leftdata->energytheta != -1);
5765  assert(rightdata->energytheta != -1);
5766  nodedata->energytheta = leftdata->energytheta + rightdata->energytheta;
5767
5768  if( leftdata->enveloplambda >= 0 )
5769  {
5770  assert(rightdata->energytheta != -1);
5771  nodedata->enveloplambda = MAX(leftdata->enveloplambda + rightdata->energytheta, rightdata->enveloplambda);
5772  }
5773  else
5774  nodedata->enveloplambda = rightdata->enveloplambda;
5775
5776  if( leftdata->enveloptheta >= 0 && rightdata->energylambda >= 0 )
5777  nodedata->enveloplambda = MAX(nodedata->enveloplambda, leftdata->enveloptheta + rightdata->energylambda);
5778
5779  SCIPdebugMsg(scip, "node <%p> lambda envelop %d\n", (void*)node, nodedata->enveloplambda);
5780
5781  if( leftdata->energylambda >= 0 && rightdata->energylambda >= 0 )
5782  {
5783  assert(rightdata->energytheta != -1);
5784  assert(leftdata->energytheta != -1);
5785  nodedata->energylambda = MAX(leftdata->energylambda + rightdata->energytheta, leftdata->energytheta + rightdata->energylambda);
5786  }
5787  else if( rightdata->energylambda >= 0 )
5788  {
5789  assert(leftdata->energytheta != -1);
5790  nodedata->energylambda = leftdata->energytheta + rightdata->energylambda;
5791  }
5792  else if( leftdata->energylambda >= 0 )
5793  {
5794  assert(rightdata->energytheta != -1);
5795  nodedata->energylambda = leftdata->energylambda + rightdata->energytheta;
5796  }
5797  else
5798  nodedata->energylambda = -1;
5799
5800  /* go to parent */
5801  node = SCIPbtnodeGetParent(node);
5802  }
5803
5804  SCIPdebugMsg(scip, "updating done\n");
5805
5806  return SCIP_OKAY;
5807 }
5808
5809 /** updates the key of the first parent on the trace which comes from left */
5810 static
5811 void updateKeyOnTrace(
5812  SCIP_BTNODE* node, /**< node to start the trace */
5813  SCIP_Real key /**< update search key */
5814  )
5815 {
5816  assert(node != NULL);
5817
5818  while( !SCIPbtnodeIsRoot(node) )
5819  {
5820  SCIP_BTNODE* parent;
5821
5822  parent = SCIPbtnodeGetParent(node);
5823  assert(parent != NULL);
5824
5825  if( SCIPbtnodeIsLeftchild(node) )
5826  {
5827  SCIP_NODEDATA* nodedata;
5828
5829  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(parent);
5830  assert(nodedata != NULL);
5831
5832  nodedata->key = key;
5833  return;
5834  }
5835
5836  node = parent;
5837  }
5838 }
5839
5840
5841 /** deletes the given node and updates all envelops */
5842 static
5844  SCIP* scip, /**< SCIP data structure */
5845  SCIP_BT* tree, /**< binary tree */
5846  SCIP_BTNODE* node /**< node to be deleted */
5847  )
5848 {
5849  SCIP_BTNODE* parent;
5850  SCIP_BTNODE* grandparent;
5851  SCIP_BTNODE* sibling;
5852
5853  assert(scip != NULL);
5854  assert(tree != NULL);
5855  assert(node != NULL);
5856
5857  assert(SCIPbtnodeIsLeaf(node));
5858  assert(!SCIPbtnodeIsRoot(node));
5859
5860  SCIPdebugMsg(scip, "delete node <%p>\n", (void*)node);
5861
5862  parent = SCIPbtnodeGetParent(node);
5863  assert(parent != NULL);
5864  if( SCIPbtnodeIsLeftchild(node) )
5865  {
5866  sibling = SCIPbtnodeGetRightchild(parent);
5867  SCIPbtnodeSetRightchild(parent, NULL);
5868  }
5869  else
5870  {
5871  sibling = SCIPbtnodeGetLeftchild(parent);
5872  SCIPbtnodeSetLeftchild(parent, NULL);
5873  }
5874  assert(sibling != NULL);
5875
5876  grandparent = SCIPbtnodeGetParent(parent);
5877
5878  if( grandparent != NULL )
5879  {
5880  /* reset parent of sibling */
5881  SCIPbtnodeSetParent(sibling, grandparent);
5882
5883  /* reset child of grandparent to sibling */
5884  if( SCIPbtnodeIsLeftchild(parent) )
5885  {
5886  SCIPbtnodeSetLeftchild(grandparent, sibling);
5887  }
5888  else
5889  {
5890  SCIP_NODEDATA* nodedata;
5891
5892  assert(SCIPbtnodeIsRightchild(parent));
5893  SCIPbtnodeSetRightchild(grandparent, sibling);
5894
5895  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(sibling);
5896
5897  updateKeyOnTrace(grandparent, nodedata->key);
5898  }
5899
5900  SCIP_CALL( updateEnvelop(scip, grandparent) );
5901  }
5902  else
5903  {
5904  SCIPbtnodeSetParent(sibling, NULL);
5905
5906  SCIPbtSetRoot(tree, sibling);
5907  }
5908
5909  SCIPbtnodeFree(tree, &parent);
5910
5911  return SCIP_OKAY;
5912 }
5913
5914 /** moves a node form the theta set into the lambda set and updates the envelops */
5915 static
5917  SCIP* scip, /**< SCIP data structure */
5918  SCIP_BT* tree, /**< binary tree */
5919  SCIP_BTNODE* node /**< node to move into the lambda set */
5920  )
5921 {
5922  SCIP_NODEDATA* nodedata;
5923
5924  assert(scip != NULL);
5925  assert(tree != NULL);
5926  assert(node != NULL);
5927
5928  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(node);
5929  assert(nodedata != NULL);
5930  assert(nodedata->intheta);
5931
5932  /* move the contributions form the theta set into the lambda set */
5933  assert(nodedata->enveloptheta != -1);
5934  assert(nodedata->energytheta != -1);
5935  assert(nodedata->enveloplambda == -1);
5936  assert(nodedata->energylambda == -1);
5937  nodedata->enveloplambda = nodedata->enveloptheta;
5938  nodedata->energylambda = nodedata->energytheta;
5939
5940  nodedata->enveloptheta = -1;
5941  nodedata->energytheta = 0;
5942  nodedata->intheta = FALSE;
5943
5944  /* update the energy and envelop values on trace */
5945  SCIP_CALL( updateEnvelop(scip, node) );
5946
5947  return SCIP_OKAY;
5948 }
5949
5950 /** inserts a node into the theta set and update the envelops */
5951 static
5953  SCIP* scip, /**< SCIP data structure */
5954  SCIP_BT* tree, /**< binary tree */
5955  SCIP_BTNODE* node, /**< node to insert */
5956  SCIP_NODEDATA** nodedatas, /**< array of node data */
5957  int* nnodedatas /**< pointer to number of node data */
5958  )
5959 {
5960  /* if the tree is empty the node will be the root node */
5961  if( SCIPbtIsEmpty(tree) )
5962  {
5963  SCIPbtSetRoot(tree, node);
5964  }
5965  else
5966  {
5967  SCIP_NODEDATA* newnodedata;
5968  SCIP_NODEDATA* leafdata;
5969  SCIP_NODEDATA* nodedata;
5970  SCIP_BTNODE* leaf;
5971  SCIP_BTNODE* newnode;
5972  SCIP_BTNODE* parent;
5973
5974  leaf = SCIPbtGetRoot(tree);
5975  assert(leaf != NULL);
5976
5977  leafdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(leaf);
5978  assert(leafdata != NULL);
5979
5980  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(node);
5981  assert(nodedata != NULL);
5982  assert(nodedata->intheta);
5983
5984  /* find the position to insert the node */
5985  while( !SCIPbtnodeIsLeaf(leaf) )
5986  {
5987  if( nodedata->key < leafdata->key )
5988  leaf = SCIPbtnodeGetLeftchild(leaf);
5989  else
5990  leaf = SCIPbtnodeGetRightchild(leaf);
5991
5992  leafdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(leaf);
5993  assert(leafdata != NULL);
5994  }
5995
5996  assert(leaf != NULL);
5997  assert(leaf != node);
5998
5999  /* create node data */
6000  SCIP_CALL( createNodedata(scip, &newnodedata) );
6001
6002  /* create a new node */
6003  SCIP_CALL( SCIPbtnodeCreate(tree, &newnode, newnodedata) );
6004  assert(newnode != NULL);
6005
6006  /* store node data to be able to delete them latter */
6007  nodedatas[*nnodedatas] = newnodedata;
6008  (*nnodedatas)++;
6009
6010  parent = SCIPbtnodeGetParent(leaf);
6011
6012  if( parent != NULL )
6013  {
6014  SCIPbtnodeSetParent(newnode, parent);
6015
6016  /* check if the node is the left child */
6017  if( SCIPbtnodeGetLeftchild(parent) == leaf )
6018  {
6019  SCIPbtnodeSetLeftchild(parent, newnode);
6020  }
6021  else
6022  {
6023  SCIPbtnodeSetRightchild(parent, newnode);
6024  }
6025  }
6026  else
6027  SCIPbtSetRoot(tree, newnode);
6028
6029  if( nodedata->key < leafdata->key )
6030  {
6031  /* node is on the left */
6032  SCIPbtnodeSetLeftchild(newnode, node);
6033  SCIPbtnodeSetRightchild(newnode, leaf);
6034  newnodedata->key = nodedata->key;
6035  }
6036  else
6037  {
6038  /* leaf is on the left */
6039  SCIPbtnodeSetLeftchild(newnode, leaf);
6040  SCIPbtnodeSetRightchild(newnode, node);
6041  newnodedata->key = leafdata->key;
6042  }
6043
6044  SCIPbtnodeSetParent(leaf, newnode);
6045  SCIPbtnodeSetParent(node, newnode);
6046  }
6047
6048  /* update envelop */
6049  SCIP_CALL( updateEnvelop(scip, node) );
6050
6051  return SCIP_OKAY;
6052 }
6053
6054 /** returns the leaf responsible for the lambda energy */
6055 static
6057  SCIP_BTNODE* node /**< node which defines the subtree beases on the lambda energy */
6058  )
6059 {
6060  SCIP_BTNODE* left;
6061  SCIP_BTNODE* right;
6062  SCIP_NODEDATA* nodedata;
6063  SCIP_NODEDATA* leftdata;
6064  SCIP_NODEDATA* rightdata;
6065
6066  assert(node != NULL);
6067
6068  nodedata = (SCIP_NODEDATA*)SCIPbtnodeGetData(node);
6069  assert(nodedata != NULL);
6070
6071  /* check if the node is the (responsible) leaf */
6072  if( SCIPbtnodeIsLeaf(node) )
6073  {
6074  assert(!nodedata->intheta);
6075  return node;
6076  }
6077
6078  left = SCIPbtnodeGetLeftchild(node);
6079  assert(left != NULL);
6080
6081  leftdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(left);
6082  assert(leftdata != NULL);
6083
6084  right = SCIPbtnodeGetRightchild(node);
6085  assert(right != NULL);
6086
6087  rightdata = (SCIP_NODEDATA*)SCIPbtnodeGetData(right);
6088  assert(rightdata != NULL);
6089
6090  assert(nodedata->energylambda != -1);
6091  assert(rightdata->energytheta != -1);
6092
6093  if( leftdata->energylambda >= 0 && nodedata->energylambda == leftdata->energylambda + rightdata->energytheta )
6095
6096  assert(leftdata->energytheta != -1);
6097  assert(rightdata->energylambda != -1);
6098  assert(nodedata->energylambda == leftdata->energytheta + rightdata->energylambda);
6099
6101 }
6102
6103 /** returns the leaf responsible for the lambda envelop */
6104 static
6106  SCIP_BTNODE* node /**< node which defines the subtree beases on the lambda envelop */
6107  )
6108 {
6109  SCIP_BTNODE* left;
6110  SCIP_BTNODE* right;
6111  SCIP_NODEDATA* nodedata;