1    	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2    	/*                                                                           */
3    	/*                  This file is part of the program and library             */
4    	/*         SCIP --- Solving Constraint Integer Programs                      */
5    	/*                                                                           */
6    	/*  Copyright (c) 2002-2023 Zuse Institute Berlin (ZIB)                      */
7    	/*                                                                           */
8    	/*  Licensed under the Apache License, Version 2.0 (the "License");          */
9    	/*  you may not use this file except in compliance with the License.         */
10   	/*  You may obtain a copy of the License at                                  */
11   	/*                                                                           */
12   	/*      http://www.apache.org/licenses/LICENSE-2.0                           */
13   	/*                                                                           */
14   	/*  Unless required by applicable law or agreed to in writing, software      */
15   	/*  distributed under the License is distributed on an "AS IS" BASIS,        */
16   	/*  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17   	/*  See the License for the specific language governing permissions and      */
18   	/*  limitations under the License.                                           */
19   	/*                                                                           */
20   	/*  You should have received a copy of the Apache-2.0 license                */
21   	/*  along with SCIP; see the file LICENSE. If not visit scipopt.org.         */
22   	/*                                                                           */
23   	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24   	
25   	/**@file   heur_rounding.c
26   	 * @ingroup DEFPLUGINS_HEUR
27   	 * @brief  LP rounding heuristic that tries to recover from intermediate infeasibilities
28   	 * @author Tobias Achterberg
29   	 */
30   	
31   	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32   	
33   	#include "blockmemshell/memory.h"
34   	#include "scip/heur_rounding.h"
35   	#include "scip/pub_heur.h"
36   	#include "scip/pub_lp.h"
37   	#include "scip/pub_message.h"
38   	#include "scip/pub_var.h"
39   	#include "scip/scip_branch.h"
40   	#include "scip/scip_heur.h"
41   	#include "scip/scip_lp.h"
42   	#include "scip/scip_mem.h"
43   	#include "scip/scip_message.h"
44   	#include "scip/scip_numerics.h"
45   	#include "scip/scip_param.h"
46   	#include "scip/scip_sol.h"
47   	#include "scip/scip_solvingstats.h"
48   	#include <string.h>
49   	
50   	#define HEUR_NAME             "rounding"
51   	#define HEUR_DESC             "LP rounding heuristic with infeasibility recovering"
52   	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_ROUNDING
53   	#define HEUR_PRIORITY         -1000
54   	#define HEUR_FREQ             1
55   	#define HEUR_FREQOFS          0
56   	#define HEUR_MAXDEPTH         -1
57   	#define HEUR_TIMING           SCIP_HEURTIMING_DURINGLPLOOP
58   	#define HEUR_USESSUBSCIP      FALSE          /**< does the heuristic use a secondary SCIP instance? */
59   	
60   	#define DEFAULT_SUCCESSFACTOR 100            /**< number of calls per found solution that are considered as standard success, 
61   	                                              * a higher factor causes the heuristic to be called more often 
62   	                                              */
63   	#define DEFAULT_ONCEPERNODE   FALSE          /**< should the heuristic only be called once per node? */
64   	
65   	/* locally defined heuristic data */
66   	struct SCIP_HeurData
67   	{
68   	   SCIP_SOL*             sol;                /**< working solution */
69   	   SCIP_Longint          lastlp;             /**< last LP number where the heuristic was applied */
70   	   int                   successfactor;      /**< number of calls per found solution that are considered as standard success, 
71   	                                              * a higher factor causes the heuristic to be called more often 
72   	                                              */
73   	   SCIP_Bool             oncepernode;        /**< should the heuristic only be called once per node? */
74   	};
75   	
76   	
77   	/*
78   	 * local methods
79   	 */
80   	
81   	/** update row violation arrays after a row's activity value changed */
82   	static
83   	void updateViolations(
84   	   SCIP*                 scip,               /**< SCIP data structure */
85   	   SCIP_ROW*             row,                /**< LP row */
86   	   SCIP_ROW**            violrows,           /**< array with currently violated rows */
87   	   int*                  violrowpos,         /**< position of LP rows in violrows array */
88   	   int*                  nviolrows,          /**< pointer to the number of currently violated rows */
89   	   SCIP_Real             oldactivity,        /**< old activity value of LP row */
90   	   SCIP_Real             newactivity         /**< new activity value of LP row */
91   	   )
92   	{
93   	   SCIP_Real lhs;
94   	   SCIP_Real rhs;
95   	   SCIP_Bool oldviol;
96   	   SCIP_Bool newviol;
97   	
98   	   assert(violrows != NULL);
99   	   assert(violrowpos != NULL);
100  	   assert(nviolrows != NULL);
101  	
102  	   lhs = SCIProwGetLhs(row);
103  	   rhs = SCIProwGetRhs(row);
104  	   oldviol = (SCIPisFeasLT(scip, oldactivity, lhs) || SCIPisFeasGT(scip, oldactivity, rhs));
105  	   newviol = (SCIPisFeasLT(scip, newactivity, lhs) || SCIPisFeasGT(scip, newactivity, rhs));
106  	   if( oldviol != newviol )
107  	   {
108  	      int rowpos;
109  	      int violpos;
110  	
111  	      rowpos = SCIProwGetLPPos(row);
112  	      assert(rowpos >= 0);
113  	
114  	      if( oldviol )
115  	      {
116  	         /* the row violation was repaired: remove row from violrows array, decrease violation count */
117  	         violpos = violrowpos[rowpos];
118  	         assert(0 <= violpos && violpos < *nviolrows);
119  	         assert(violrows[violpos] == row);
120  	         violrowpos[rowpos] = -1;
121  	         if( violpos != *nviolrows-1 )
122  	         {
123  	            violrows[violpos] = violrows[*nviolrows-1];
124  	            violrowpos[SCIProwGetLPPos(violrows[violpos])] = violpos;
125  	         }
126  	         (*nviolrows)--;
127  	      }
128  	      else
129  	      {
130  	         /* the row is now violated: add row to violrows array, increase violation count */
131  	         assert(violrowpos[rowpos] == -1);
132  	         violrows[*nviolrows] = row;
133  	         violrowpos[rowpos] = *nviolrows;
134  	         (*nviolrows)++;
135  	      }
136  	   }
137  	}
138  	
139  	/** update row activities after a variable's solution value changed */
140  	static
141  	SCIP_RETCODE updateActivities(
142  	   SCIP*                 scip,               /**< SCIP data structure */
143  	   SCIP_Real*            activities,         /**< LP row activities */
144  	   SCIP_ROW**            violrows,           /**< array with currently violated rows */
145  	   int*                  violrowpos,         /**< position of LP rows in violrows array */
146  	   int*                  nviolrows,          /**< pointer to the number of currently violated rows */
147  	   int                   nlprows,            /**< number of rows in current LP */
148  	   SCIP_VAR*             var,                /**< variable that has been changed */
149  	   SCIP_Real             oldsolval,          /**< old solution value of variable */
150  	   SCIP_Real             newsolval           /**< new solution value of variable */
151  	   )
152  	{
153  	   SCIP_COL* col;
154  	   SCIP_ROW** colrows;
155  	   SCIP_Real* colvals;
156  	   SCIP_Real delta;
157  	   int ncolrows;
158  	   int r;
159  	
160  	   assert(activities != NULL);
161  	   assert(nviolrows != NULL);
162  	   assert(0 <= *nviolrows && *nviolrows <= nlprows);
163  	
164  	   delta = newsolval - oldsolval;
165  	   col = SCIPvarGetCol(var);
166  	   colrows = SCIPcolGetRows(col);
167  	   colvals = SCIPcolGetVals(col);
168  	   ncolrows = SCIPcolGetNLPNonz(col);
169  	   assert(ncolrows == 0 || (colrows != NULL && colvals != NULL));
170  	
171  	   for( r = 0; r < ncolrows; ++r )
172  	   {
173  	      SCIP_ROW* row;
174  	      int rowpos;
175  	
176  	      row = colrows[r];
177  	      rowpos = SCIProwGetLPPos(row);
178  	      assert(-1 <= rowpos && rowpos < nlprows);
179  	
180  	      if( rowpos >= 0 && !SCIProwIsLocal(row) )
181  	      {
182  	         SCIP_Real oldactivity;
183  	         SCIP_Real newactivity;
184  	
185  	         assert(SCIProwIsInLP(row));
186  	
187  	         /* update row activity */
188  	         oldactivity = activities[rowpos];
189  	         if( !SCIPisInfinity(scip, -oldactivity) && !SCIPisInfinity(scip, oldactivity) )
190  	         {
191  	            newactivity = oldactivity + delta * colvals[r];
192  	            if( SCIPisInfinity(scip, newactivity) )
193  	               newactivity = SCIPinfinity(scip);
194  	            else if( SCIPisInfinity(scip, -newactivity) )
195  	               newactivity = -SCIPinfinity(scip);
196  	            activities[rowpos] = newactivity;
197  	
198  	            /* update row violation arrays */
199  	            updateViolations(scip, row, violrows, violrowpos, nviolrows, oldactivity, newactivity);
200  	         }
201  	      }
202  	   }
203  	
204  	   return SCIP_OKAY;
205  	}
206  	
207  	/** returns a variable, that pushes activity of the row in the given direction with minimal negative impact on other rows;
208  	 *  if variables have equal impact, chooses the one with best objective value improvement in corresponding direction;
209  	 *  rounding in a direction is forbidden, if this forces the objective value over the upper bound
210  	 */
211  	static
212  	SCIP_RETCODE selectRounding(
213  	   SCIP*                 scip,               /**< SCIP data structure */
214  	   SCIP_SOL*             sol,                /**< primal solution */
215  	   SCIP_Real             minobj,             /**< minimal objective value possible after rounding remaining fractional vars */
216  	   SCIP_ROW*             row,                /**< LP row */
217  	   int                   direction,          /**< should the activity be increased (+1) or decreased (-1)? */
218  	   SCIP_VAR**            roundvar,           /**< pointer to store the rounding variable, returns NULL if impossible */
219  	   SCIP_Real*            oldsolval,          /**< pointer to store old (fractional) solution value of rounding variable */
220  	   SCIP_Real*            newsolval           /**< pointer to store new (rounded) solution value of rounding variable */
221  	   )
222  	{
223  	   SCIP_COL* col;
224  	   SCIP_VAR* var;
225  	   SCIP_Real val;
226  	   SCIP_COL** rowcols;
227  	   SCIP_Real* rowvals;
228  	   SCIP_Real solval;
229  	   SCIP_Real roundval;
230  	   SCIP_Real obj;
231  	   SCIP_Real deltaobj;
232  	   SCIP_Real bestdeltaobj;
233  	   SCIP_VARTYPE vartype;
234  	   int nrowcols;
235  	   int nlocks;
236  	   int minnlocks;
237  	   int c;
238  	
239  	   assert(direction == +1 || direction == -1);
240  	   assert(roundvar != NULL);
241  	   assert(oldsolval != NULL);
242  	   assert(newsolval != NULL);
243  	
244  	   /* get row entries */
245  	   rowcols = SCIProwGetCols(row);
246  	   rowvals = SCIProwGetVals(row);
247  	   nrowcols = SCIProwGetNLPNonz(row);
248  	
249  	   /* select rounding variable */
250  	   minnlocks = INT_MAX;
251  	   bestdeltaobj = SCIPinfinity(scip);
252  	   *roundvar = NULL;
253  	   for( c = 0; c < nrowcols; ++c )
254  	   {
255  	      col = rowcols[c];
256  	      var = SCIPcolGetVar(col);
257  	
258  	      vartype = SCIPvarGetType(var);
259  	      if( vartype == SCIP_VARTYPE_BINARY || vartype == SCIP_VARTYPE_INTEGER )
260  	      {
261  	         solval = SCIPgetSolVal(scip, sol, var);
262  	
263  	         if( !SCIPisFeasIntegral(scip, solval) )
264  	         {
265  	            val = rowvals[c];
266  	            obj = SCIPvarGetObj(var);
267  	
268  	            if( direction * val < 0.0 )
269  	            {
270  	               /* rounding down */
271  	               nlocks = SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL);
272  	               if( nlocks <= minnlocks )
273  	               {
274  	                  roundval = SCIPfeasFloor(scip, solval);
275  	                  deltaobj = obj * (roundval - solval);
276  	                  if( (nlocks < minnlocks || deltaobj < bestdeltaobj) && minobj - obj < SCIPgetCutoffbound(scip) )
277  	                  {
278  	                     minnlocks = nlocks;
279  	                     bestdeltaobj = deltaobj;
280  	                     *roundvar = var;
281  	                     *oldsolval = solval;
282  	                     *newsolval = roundval;
283  	                  }
284  	               }
285  	            }
286  	            else
287  	            {
288  	               /* rounding up */
289  	               assert(direction * val > 0.0);
290  	               nlocks = SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL);
291  	               if( nlocks <= minnlocks )
292  	               {
293  	                  roundval = SCIPfeasCeil(scip, solval);
294  	                  deltaobj = obj * (roundval - solval);
295  	                  if( (nlocks < minnlocks || deltaobj < bestdeltaobj) && minobj + obj < SCIPgetCutoffbound(scip) )
296  	                  {
297  	                     minnlocks = nlocks;
298  	                     bestdeltaobj = deltaobj;
299  	                     *roundvar = var;
300  	                     *oldsolval = solval;
301  	                     *newsolval = roundval;
302  	                  }
303  	               }
304  	            }
305  	         }
306  	      }
307  	   }
308  	
309  	   return SCIP_OKAY;
310  	}
311  	
312  	/** returns a variable, that increases the activity of the row */
313  	static
314  	SCIP_RETCODE selectIncreaseRounding(
315  	   SCIP*                 scip,               /**< SCIP data structure */
316  	   SCIP_SOL*             sol,                /**< primal solution */
317  	   SCIP_Real             minobj,             /**< minimal objective value possible after rounding remaining fractional vars */
318  	   SCIP_ROW*             row,                /**< LP row */
319  	   SCIP_VAR**            roundvar,           /**< pointer to store the rounding variable, returns NULL if impossible */
320  	   SCIP_Real*            oldsolval,          /**< old (fractional) solution value of rounding variable */
321  	   SCIP_Real*            newsolval           /**< new (rounded) solution value of rounding variable */
322  	   )
323  	{
324  	   return selectRounding(scip, sol, minobj, row, +1, roundvar, oldsolval, newsolval);
325  	}
326  	
327  	/** returns a variable, that decreases the activity of the row */
328  	static
329  	SCIP_RETCODE selectDecreaseRounding(
330  	   SCIP*                 scip,               /**< SCIP data structure */
331  	   SCIP_SOL*             sol,                /**< primal solution */
332  	   SCIP_Real             minobj,             /**< minimal objective value possible after rounding remaining fractional vars */
333  	   SCIP_ROW*             row,                /**< LP row */
334  	   SCIP_VAR**            roundvar,           /**< pointer to store the rounding variable, returns NULL if impossible */
335  	   SCIP_Real*            oldsolval,          /**< old (fractional) solution value of rounding variable */
336  	   SCIP_Real*            newsolval           /**< new (rounded) solution value of rounding variable */
337  	   )
338  	{
339  	   return selectRounding(scip, sol, minobj, row, -1, roundvar, oldsolval, newsolval);
340  	}
341  	
342  	/** returns a fractional variable, that has most impact on rows in opposite direction, i.e. that is most crucial to
343  	 *  fix in the other direction;
344  	 *  if variables have equal impact, chooses the one with best objective value improvement in corresponding direction;
345  	 *  rounding in a direction is forbidden, if this forces the objective value over the upper bound
346  	 */
347  	static
348  	SCIP_RETCODE selectEssentialRounding(
349  	   SCIP*                 scip,               /**< SCIP data structure */
350  	   SCIP_SOL*             sol,                /**< primal solution */
351  	   SCIP_Real             minobj,             /**< minimal objective value possible after rounding remaining fractional vars */
352  	   SCIP_VAR**            lpcands,            /**< fractional variables in LP */
353  	   int                   nlpcands,           /**< number of fractional variables in LP */
354  	   SCIP_VAR**            roundvar,           /**< pointer to store the rounding variable, returns NULL if impossible */
355  	   SCIP_Real*            oldsolval,          /**< old (fractional) solution value of rounding variable */
356  	   SCIP_Real*            newsolval           /**< new (rounded) solution value of rounding variable */
357  	   )
358  	{
359  	   SCIP_VAR* var;
360  	   SCIP_Real solval;
361  	   SCIP_Real roundval;
362  	   SCIP_Real obj;
363  	   SCIP_Real deltaobj;
364  	   SCIP_Real bestdeltaobj;
365  	   int maxnlocks;
366  	   int nlocks;
367  	   int v;
368  	
369  	   assert(roundvar != NULL);
370  	   assert(oldsolval != NULL);
371  	   assert(newsolval != NULL);
372  	
373  	   /* select rounding variable */
374  	   maxnlocks = -1;
375  	   bestdeltaobj = SCIPinfinity(scip);
376  	   *roundvar = NULL;
377  	   for( v = 0; v < nlpcands; ++v )
378  	   {
379  	      var = lpcands[v];
380  	      assert(SCIPvarGetType(var) == SCIP_VARTYPE_BINARY || SCIPvarGetType(var) == SCIP_VARTYPE_INTEGER);
381  	
382  	      solval = SCIPgetSolVal(scip, sol, var);
383  	      if( !SCIPisFeasIntegral(scip, solval) )
384  	      {
385  	         obj = SCIPvarGetObj(var);
386  	
387  	         /* rounding down */
388  	         nlocks = SCIPvarGetNLocksUpType(var, SCIP_LOCKTYPE_MODEL);
389  	         if( nlocks >= maxnlocks )
390  	         {
391  	            roundval = SCIPfeasFloor(scip, solval);
392  	            deltaobj = obj * (roundval - solval);
393  	            if( (nlocks > maxnlocks || deltaobj < bestdeltaobj) && minobj - obj < SCIPgetCutoffbound(scip) )
394  	            {
395  	               maxnlocks = nlocks;
396  	               bestdeltaobj = deltaobj;
397  	               *roundvar = var;
398  	               *oldsolval = solval;
399  	               *newsolval = roundval;
400  	            }
401  	         }
402  	
403  	         /* rounding up */
404  	         nlocks = SCIPvarGetNLocksDownType(var, SCIP_LOCKTYPE_MODEL);
405  	         if( nlocks >= maxnlocks )
406  	         {
407  	            roundval = SCIPfeasCeil(scip, solval);
408  	            deltaobj = obj * (roundval - solval);
409  	            if( (nlocks > maxnlocks || deltaobj < bestdeltaobj) && minobj + obj < SCIPgetCutoffbound(scip) )
410  	            {
411  	               maxnlocks = nlocks;
412  	               bestdeltaobj = deltaobj;
413  	               *roundvar = var;
414  	               *oldsolval = solval;
415  	               *newsolval = roundval;
416  	            }
417  	         }
418  	      }
419  	   }
420  	
421  	   return SCIP_OKAY;
422  	}
423  	
424  	
425  	/*
426  	 * Callback methods
427  	 */
428  	
429  	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
430  	static
431  	SCIP_DECL_HEURCOPY(heurCopyRounding)
432  	{  /*lint --e{715}*/
433  	   assert(scip != NULL);
434  	   assert(heur != NULL);
435  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
436  	
437  	   /* call inclusion method of primal heuristic */
438  	   SCIP_CALL( SCIPincludeHeurRounding(scip) );
439  	
440  	   return SCIP_OKAY;
441  	}
442  	
443  	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
444  	static
445  	SCIP_DECL_HEURFREE(heurFreeRounding) /*lint --e{715}*/
446  	{  /*lint --e{715}*/
447  	   SCIP_HEURDATA* heurdata;
448  	
449  	   assert(heur != NULL);
450  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
451  	   assert(scip != NULL);
452  	
453  	   /* free heuristic data */
454  	   heurdata = SCIPheurGetData(heur);
455  	   assert(heurdata != NULL);
456  	   SCIPfreeBlockMemory(scip, &heurdata);
457  	   SCIPheurSetData(heur, NULL);
458  	
459  	   return SCIP_OKAY;
460  	}
461  	
462  	/** initialization method of primal heuristic (called after problem was transformed) */
463  	static
464  	SCIP_DECL_HEURINIT(heurInitRounding) /*lint --e{715}*/
465  	{  /*lint --e{715}*/
466  	   SCIP_HEURDATA* heurdata;
467  	
468  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
469  	
470  	   heurdata = SCIPheurGetData(heur);
471  	   assert(heurdata != NULL);
472  	
473  	   /* create heuristic data */
474  	   SCIP_CALL( SCIPcreateSol(scip, &heurdata->sol, heur) );
475  	   heurdata->lastlp = -1;
476  	
477  	   return SCIP_OKAY;
478  	}
479  	
480  	/** deinitialization method of primal heuristic (called before transformed problem is freed) */
481  	static
482  	SCIP_DECL_HEUREXIT(heurExitRounding) /*lint --e{715}*/
483  	{  /*lint --e{715}*/
484  	   SCIP_HEURDATA* heurdata;
485  	
486  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
487  	
488  	   /* free heuristic data */
489  	   heurdata = SCIPheurGetData(heur);
490  	   assert(heurdata != NULL);
491  	   SCIP_CALL( SCIPfreeSol(scip, &heurdata->sol) );
492  	
493  	   return SCIP_OKAY;
494  	}
495  	
496  	/** solving process initialization method of primal heuristic (called when branch and bound process is about to begin) */
497  	static
498  	SCIP_DECL_HEURINITSOL(heurInitsolRounding)
499  	{
500  	   SCIP_HEURDATA* heurdata;
501  	
502  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
503  	
504  	   heurdata = SCIPheurGetData(heur);
505  	   assert(heurdata != NULL);
506  	   heurdata->lastlp = -1;
507  	
508  	   /* change the heuristic's timingmask, if nit should be called only once per node */
509  	   if( heurdata->oncepernode )
510  	      SCIPheurSetTimingmask(heur, SCIP_HEURTIMING_AFTERLPNODE);
511  	
512  	   return SCIP_OKAY;
513  	}
514  	
515  	/** solving process deinitialization method of primal heuristic (called before branch and bound process data is freed) */
516  	static
517  	SCIP_DECL_HEUREXITSOL(heurExitsolRounding)
518  	{
519  	   /* reset the timing mask to its default value */
520  	   SCIPheurSetTimingmask(heur, HEUR_TIMING);
521  	
522  	   return SCIP_OKAY;
523  	}
524  	
525  	
526  	/** execution method of primal heuristic */
527  	static
528  	SCIP_DECL_HEUREXEC(heurExecRounding) /*lint --e{715}*/
529  	{  /*lint --e{715}*/
530  	   SCIP_HEURDATA* heurdata;
531  	   SCIP_SOL* sol;
532  	   SCIP_VAR** lpcands;
533  	   SCIP_Real* lpcandssol;
534  	   SCIP_ROW** lprows;
535  	   SCIP_Real* activities;
536  	   SCIP_ROW** violrows;
537  	   int* violrowpos;
538  	   SCIP_Real obj;
539  	   SCIP_Real bestroundval;
540  	   SCIP_Real minobj;
541  	   int nlpcands;
542  	   int nlprows;
543  	   int nfrac;
544  	   int nviolrows;
545  	   int c;
546  	   int r;
547  	   SCIP_Longint nlps;
548  	   SCIP_Longint ncalls;
549  	   SCIP_Longint nsolsfound;
550  	   SCIP_Longint nnodes;
551  	
552  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
553  	   assert(scip != NULL);
554  	   assert(result != NULL);
555  	   assert(SCIPhasCurrentNodeLP(scip));
556  	
557  	   *result = SCIP_DIDNOTRUN;
558  	
559  	   /* only call heuristic, if an optimal LP solution is at hand */
560  	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
561  	      return SCIP_OKAY;
562  	
563  	   /* only call heuristic, if the LP objective value is smaller than the cutoff bound */
564  	   if( SCIPisGE(scip, SCIPgetLPObjval(scip), SCIPgetCutoffbound(scip)) )
565  	      return SCIP_OKAY;
566  	
567  	   /* get heuristic data */
568  	   heurdata = SCIPheurGetData(heur);
569  	   assert(heurdata != NULL);
570  	
571  	   /* don't call heuristic, if we have already processed the current LP solution */
572  	   nlps = SCIPgetNLPs(scip);
573  	   if( nlps == heurdata->lastlp )
574  	      return SCIP_OKAY;
575  	   heurdata->lastlp = nlps;
576  	
577  	   /* don't call heuristic, if it was not successful enough in the past */
578  	   ncalls = SCIPheurGetNCalls(heur);
579  	   nsolsfound = 10*SCIPheurGetNBestSolsFound(heur) + SCIPheurGetNSolsFound(heur);
580  	   nnodes = SCIPgetNNodes(scip);
581  	   if( nnodes % ((ncalls/heurdata->successfactor)/(nsolsfound+1)+1) != 0 )
582  	      return SCIP_OKAY;
583  	
584  	   /* get fractional variables, that should be integral */
585  	   SCIP_CALL( SCIPgetLPBranchCands(scip, &lpcands, &lpcandssol, NULL, &nlpcands, NULL, NULL) );
586  	   nfrac = nlpcands;
587  	
588  	   /* only call heuristic, if LP solution is fractional */
589  	   if( nfrac == 0 )
590  	      return SCIP_OKAY;
591  	
592  	   *result = SCIP_DIDNOTFIND;
593  	
594  	   /* get LP rows */
595  	   SCIP_CALL( SCIPgetLPRowsData(scip, &lprows, &nlprows) );
596  	
597  	   SCIPdebugMsg(scip, "executing rounding heuristic: %d LP rows, %d fractionals\n", nlprows, nfrac);
598  	
599  	   /* get memory for activities, violated rows, and row violation positions */
600  	   SCIP_CALL( SCIPallocBufferArray(scip, &activities, nlprows) );
601  	   SCIP_CALL( SCIPallocBufferArray(scip, &violrows, nlprows) );
602  	   SCIP_CALL( SCIPallocBufferArray(scip, &violrowpos, nlprows) );
603  	
604  	   /* get the activities for all globally valid rows;
605  	    * the rows should be feasible, but due to numerical inaccuracies in the LP solver, they can be violated
606  	    */
607  	   nviolrows = 0;
608  	   for( r = 0; r < nlprows; ++r )
609  	   {
610  	      SCIP_ROW* row;
611  	
612  	      row = lprows[r];
613  	      assert(SCIProwGetLPPos(row) == r);
614  	
615  	      if( !SCIProwIsLocal(row) )
616  	      {
617  	         activities[r] = SCIPgetRowActivity(scip, row);
618  	         if( SCIPisFeasLT(scip, activities[r], SCIProwGetLhs(row))
619  	            || SCIPisFeasGT(scip, activities[r], SCIProwGetRhs(row)) )
620  	         {
621  	            violrows[nviolrows] = row;
622  	            violrowpos[r] = nviolrows;
623  	            nviolrows++;
624  	         }
625  	         else
626  	            violrowpos[r] = -1;
627  	      }
628  	   }
629  	
630  	   /* get the working solution from heuristic's local data */
631  	   sol = heurdata->sol;
632  	   assert(sol != NULL);
633  	
634  	   /* copy the current LP solution to the working solution */
635  	   SCIP_CALL( SCIPlinkLPSol(scip, sol) );
636  	
637  	   /* calculate the minimal objective value possible after rounding fractional variables */
638  	   minobj = SCIPgetSolTransObj(scip, sol);
639  	   assert(minobj < SCIPgetCutoffbound(scip));
640  	   for( c = 0; c < nlpcands; ++c )
641  	   {
642  	      obj = SCIPvarGetObj(lpcands[c]);
643  	      bestroundval = obj > 0.0 ? SCIPfeasFloor(scip, lpcandssol[c]) : SCIPfeasCeil(scip, lpcandssol[c]);
644  	      minobj += obj * (bestroundval - lpcandssol[c]);
645  	   }
646  	
647  	   /* try to round remaining variables in order to become/stay feasible */
648  	   while( nfrac > 0 )
649  	   {
650  	      SCIP_VAR* roundvar;
651  	      SCIP_Real oldsolval;
652  	      SCIP_Real newsolval;
653  	
654  	      SCIPdebugMsg(scip, "rounding heuristic: nfrac=%d, nviolrows=%d, obj=%g (best possible obj: %g)\n",
655  	         nfrac, nviolrows, SCIPgetSolOrigObj(scip, sol), SCIPretransformObj(scip, minobj));
656  	
657  	      /* minobj < SCIPgetCutoffbound(scip) should be true, otherwise the rounding variable selection
658  	       * should have returned NULL. Due to possible cancellation we use SCIPisLE. */
659  	      assert( SCIPisLE(scip, minobj, SCIPgetCutoffbound(scip)) );
660  	
661  	      /* choose next variable to process:
662  	       *  - if a violated row exists, round a variable decreasing the violation, that has least impact on other rows
663  	       *  - otherwise, round a variable, that has strongest devastating impact on rows in opposite direction
664  	       */
665  	      if( nviolrows > 0 )
666  	      {
667  	         SCIP_ROW* row;
668  	         int rowpos;
669  	
670  	         row = violrows[nviolrows-1];
671  	         rowpos = SCIProwGetLPPos(row);
672  	         assert(0 <= rowpos && rowpos < nlprows);
673  	         assert(violrowpos[rowpos] == nviolrows-1);
674  	
675  	         SCIPdebugMsg(scip, "rounding heuristic: try to fix violated row <%s>: %g <= %g <= %g\n",
676  	            SCIProwGetName(row), SCIProwGetLhs(row), activities[rowpos], SCIProwGetRhs(row));
677  	         if( SCIPisFeasLT(scip, activities[rowpos], SCIProwGetLhs(row)) )
678  	         {
679  	            /* lhs is violated: select a variable rounding, that increases the activity */
680  	            SCIP_CALL( selectIncreaseRounding(scip, sol, minobj, row, &roundvar, &oldsolval, &newsolval) );
681  	         }
682  	         else
683  	         {
684  	            assert(SCIPisFeasGT(scip, activities[rowpos], SCIProwGetRhs(row)));
685  	            /* rhs is violated: select a variable rounding, that decreases the activity */
686  	            SCIP_CALL( selectDecreaseRounding(scip, sol, minobj, row, &roundvar, &oldsolval, &newsolval) );
687  	         }
688  	      }
689  	      else
690  	      {
691  	         SCIPdebugMsg(scip, "rounding heuristic: search rounding variable and try to stay feasible\n");
692  	         SCIP_CALL( selectEssentialRounding(scip, sol, minobj, lpcands, nlpcands, &roundvar, &oldsolval, &newsolval) );
693  	      }
694  	
695  	      /* check, whether rounding was possible */
696  	      if( roundvar == NULL )
697  	      {
698  	         SCIPdebugMsg(scip, "rounding heuristic:  -> didn't find a rounding variable\n");
699  	         break;
700  	      }
701  	
702  	      SCIPdebugMsg(scip, "rounding heuristic:  -> round var <%s>, oldval=%g, newval=%g, obj=%g\n",
703  	         SCIPvarGetName(roundvar), oldsolval, newsolval, SCIPvarGetObj(roundvar));
704  	
705  	      /* update row activities of globally valid rows */
706  	      SCIP_CALL( updateActivities(scip, activities, violrows, violrowpos, &nviolrows, nlprows, 
707  	            roundvar, oldsolval, newsolval) );
708  	
709  	      /* store new solution value and decrease fractionality counter */
710  	      SCIP_CALL( SCIPsetSolVal(scip, sol, roundvar, newsolval) );
711  	      nfrac--;
712  	
713  	      /* update minimal objective value possible after rounding remaining variables */
714  	      obj = SCIPvarGetObj(roundvar);
715  	      if( obj > 0.0 && newsolval > oldsolval )
716  	         minobj += obj;
717  	      else if( obj < 0.0 && newsolval < oldsolval )
718  	         minobj -= obj;
719  	
720  	      SCIPdebugMsg(scip, "rounding heuristic:  -> nfrac=%d, nviolrows=%d, obj=%g (best possible obj: %g)\n",
721  	         nfrac, nviolrows, SCIPgetSolOrigObj(scip, sol), SCIPretransformObj(scip, minobj));
722  	   }
723  	
724  	   /* check, if the new solution is feasible */
725  	   if( nfrac == 0 && nviolrows == 0 )
726  	   {
727  	      SCIP_Bool stored;
728  	
729  	      /* check solution for feasibility, and add it to solution store if possible
730  	       * neither integrality nor feasibility of LP rows has to be checked, because this is already
731  	       * done in the rounding heuristic itself; however, be better check feasibility of LP rows,
732  	       * because of numerical problems with activity updating
733  	       */
734  	      SCIP_CALL( SCIPtrySol(scip, sol, FALSE, FALSE, FALSE, FALSE, TRUE, &stored) );
735  	
736  	      if( stored )
737  	      {
738  	#ifdef SCIP_DEBUG
739  	         SCIPdebugMsg(scip, "found feasible rounded solution:\n");
740  	         SCIP_CALL( SCIPprintSol(scip, sol, NULL, FALSE) );
741  	#endif
742  	         *result = SCIP_FOUNDSOL;
743  	      }
744  	   }
745  	
746  	   /* free memory buffers */
747  	   SCIPfreeBufferArray(scip, &violrowpos);
748  	   SCIPfreeBufferArray(scip, &violrows);
749  	   SCIPfreeBufferArray(scip, &activities);
750  	
751  	   return SCIP_OKAY;
752  	}
753  	
754  	
755  	/*
756  	 * heuristic specific interface methods
757  	 */
758  	
759  	/** creates the rounding heuristic with infeasibility recovering and includes it in SCIP */
760  	SCIP_RETCODE SCIPincludeHeurRounding(
761  	   SCIP*                 scip                /**< SCIP data structure */
762  	   )
763  	{
764  	   SCIP_HEURDATA* heurdata;
765  	   SCIP_HEUR* heur;
766  	
767  	   /* create Rounding primal heuristic data */
768  	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
769  	
770  	   /* include primal heuristic */
771  	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
772  	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
773  	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecRounding, heurdata) );
774  	
775  	   assert(heur != NULL);
776  	
777  	   /* set non-NULL pointers to callback methods */
778  	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyRounding) );
779  	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeRounding) );
780  	   SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitRounding) );
781  	   SCIP_CALL( SCIPsetHeurExit(scip, heur, heurExitRounding) );
782  	   SCIP_CALL( SCIPsetHeurInitsol(scip, heur, heurInitsolRounding) );
783  	   SCIP_CALL( SCIPsetHeurExitsol(scip, heur, heurExitsolRounding) );
784  	
785  	   /* add rounding primal heuristic parameters */
786  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/successfactor",
787  	         "number of calls per found solution that are considered as standard success, a higher factor causes the heuristic to be called more often",
788  	         &heurdata->successfactor, TRUE, DEFAULT_SUCCESSFACTOR, -1, INT_MAX, NULL, NULL) );
789  	
790  	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/oncepernode",
791  	         "should the heuristic only be called once per node?",
792  	         &heurdata->oncepernode, TRUE, DEFAULT_ONCEPERNODE, NULL, NULL) );
793  	
794  	   return SCIP_OKAY;
795  	}
796