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_distributiondiving.c
26   	 * @ingroup DEFPLUGINS_HEUR
27   	 * @brief Diving heuristic that chooses fixings w.r.t. changes in the solution density after Pryor and Chinneck.
28   	 * @author Gregor Hendel
29   	 */
30   	
31   	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32   	
33   	#include "blockmemshell/memory.h"
34   	#include "scip/branch_distribution.h"
35   	#include "scip/heur_distributiondiving.h"
36   	#include "scip/heuristics.h"
37   	#include "scip/pub_event.h"
38   	#include "scip/pub_heur.h"
39   	#include "scip/pub_lp.h"
40   	#include "scip/pub_message.h"
41   	#include "scip/pub_var.h"
42   	#include "scip/scip_event.h"
43   	#include "scip/scip_general.h"
44   	#include "scip/scip_heur.h"
45   	#include "scip/scip_lp.h"
46   	#include "scip/scip_mem.h"
47   	#include "scip/scip_message.h"
48   	#include "scip/scip_numerics.h"
49   	#include "scip/scip_param.h"
50   	#include "scip/scip_prob.h"
51   	#include "scip/scip_probing.h"
52   	#include "scip/scip_sol.h"
53   	#include <string.h>
54   	
55   	#define HEUR_NAME             "distributiondiving"
56   	#define HEUR_DESC             "Diving heuristic that chooses fixings w.r.t. changes in the solution density"
57   	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_DIVING
58   	#define HEUR_PRIORITY         -1003300
59   	#define HEUR_FREQ             10
60   	#define HEUR_FREQOFS          3
61   	#define HEUR_MAXDEPTH         -1
62   	#define HEUR_TIMING           SCIP_HEURTIMING_AFTERLPPLUNGE
63   	#define HEUR_USESSUBSCIP      FALSE  /**< does the heuristic use a secondary SCIP instance? */
64   	#define EVENT_DISTRIBUTION    SCIP_EVENTTYPE_BOUNDCHANGED /**< the event type to be handled by this event handler */
65   	#define EVENTHDLR_NAME        "eventhdlr_distributiondiving"
66   	#define DIVESET_DIVETYPES     SCIP_DIVETYPE_INTEGRALITY /**< bit mask that represents all supported dive types */
67   	#define DIVESET_ISPUBLIC      FALSE  /**< is this dive set publicly available (ie., can be used by other primal heuristics?) */
68   	
69   	#define SQUARED(x) ((x) * (x))
70   	/*
71   	 * Default parameter settings
72   	 */
73   	
74   	#define DEFAULT_MINRELDEPTH         0.0 /**< minimal relative depth to start diving */
75   	#define DEFAULT_MAXRELDEPTH         1.0 /**< maximal relative depth to start diving */
76   	#define DEFAULT_MAXLPITERQUOT      0.05 /**< maximal fraction of diving LP iterations compared to node LP iterations */
77   	#define DEFAULT_MAXLPITEROFS       1000 /**< additional number of allowed LP iterations */
78   	#define DEFAULT_MAXDIVEUBQUOT       0.8 /**< maximal quotient (curlowerbound - lowerbound)/(cutoffbound - lowerbound)
79   	                                         *   where diving is performed (0.0: no limit) */
80   	#define DEFAULT_MAXDIVEAVGQUOT      0.0 /**< maximal quotient (curlowerbound - lowerbound)/(avglowerbound - lowerbound)
81   	                                         *   where diving is performed (0.0: no limit) */
82   	#define DEFAULT_MAXDIVEUBQUOTNOSOL  0.1 /**< maximal UBQUOT when no solution was found yet (0.0: no limit) */
83   	#define DEFAULT_MAXDIVEAVGQUOTNOSOL 0.0 /**< maximal AVGQUOT when no solution was found yet (0.0: no limit) */
84   	#define DEFAULT_BACKTRACK          TRUE /**< use one level of backtracking if infeasibility is encountered? */
85   	#define DEFAULT_LPRESOLVEDOMCHGQUOT 0.15 /**< percentage of immediate domain changes during probing to trigger LP resolve */
86   	#define DEFAULT_LPSOLVEFREQ           0 /**< LP solve frequency for diving heuristics */
87   	#define DEFAULT_ONLYLPBRANCHCANDS  TRUE /**< should only LP branching candidates be considered instead of the slower but
88   	                                         *   more general constraint handler diving variable selection? */
89   	
90   	#define SCOREPARAM_VALUES "lhwvd"       /**< the score;largest 'd'ifference, 'l'owest cumulative probability,'h'ighest c.p.,
91   	                                         *   'v'otes lowest c.p., votes highest c.p.('w'), 'r'evolving */
92   	#define SCOREPARAM_VALUESLEN 5
93   	#define DEFAULT_SCOREPARAM 'r'          /**< default scoring parameter to guide the diving */
94   	#define DEFAULT_RANDSEED   117          /**< initial seed for random number generation */
95   	
96   	/* locally defined heuristic data */
97   	struct SCIP_HeurData
98   	{
99   	   SCIP_SOL*             sol;                /**< working solution */
100  	   SCIP_EVENTHDLR*       eventhdlr;          /**< event handler pointer */
101  	   SCIP_VAR**            updatedvars;        /**< variables to process bound change events for */
102  	   SCIP_Real*            rowmeans;           /**< row activity mean values for all rows */
103  	   SCIP_Real*            rowvariances;       /**< row activity variances for all rows */
104  	   SCIP_Real*            currentubs;         /**< variable upper bounds as currently saved in the row activities */
105  	   SCIP_Real*            currentlbs;         /**< variable lower bounds as currently saved in the row activities */
106  	   int*                  rowinfinitiesdown;  /**< count the number of variables with infinite bounds which allow for always
107  	                                              *   repairing the constraint right hand side */
108  	   int*                  rowinfinitiesup;    /**< count the number of variables with infinite bounds which allow for always
109  	                                              *   repairing the constraint left hand side */
110  	   int*                  varposs;            /**< array of variable positions in the updated variables array */
111  	   int*                  varfilterposs;      /**< array of event filter positions for variable events */
112  	   int                   nupdatedvars;       /**< the current number of variables with pending bound changes */
113  	   int                   memsize;            /**< memory size of current arrays, needed for dynamic reallocation */
114  	   int                   varpossmemsize;     /**< memory size of updated vars and varposs array */
115  	
116  	   char                  scoreparam;         /**< score user parameter */
117  	   char                  score;              /**< score to be used depending on user parameter to use fixed score or revolve */
118  	};
119  	
120  	struct SCIP_EventhdlrData
121  	{
122  	   SCIP_HEURDATA*  heurdata;     /**< the heuristic data to access distribution arrays */
123  	};
124  	/*
125  	 * local methods
126  	 */
127  	
128  	/** ensure that maxindex + 1 rows can be represented in data arrays; memory gets reallocated with 10% extra space
129  	 *  to save some time for future allocations */
130  	static
131  	SCIP_RETCODE heurdataEnsureArraySize(
132  	   SCIP*                 scip,               /**< SCIP data structure */
133  	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
134  	   int                   maxindex            /**< row index at hand (size must be at least this large) */
135  	   )
136  	{
137  	   int newsize;
138  	   int r;
139  	
140  	   /* maxindex fits in current array -> nothing to do */
141  	   if( maxindex < heurdata->memsize )
142  	      return SCIP_OKAY;
143  	
144  	   /* new memory size is the max index + 1 plus 10% additional space */
145  	   newsize = (int)SCIPfeasCeil(scip, (maxindex + 1) * 1.1);
146  	   assert(newsize > heurdata->memsize);
147  	   assert(heurdata->memsize >= 0);
148  	
149  	   /* alloc memory arrays for row information */
150  	   if( heurdata->memsize == 0 )
151  	   {
152  	      SCIP_VAR** vars;
153  	      int v;
154  	      int nvars;
155  	
156  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->rowinfinitiesdown, newsize) );
157  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->rowinfinitiesup, newsize) );
158  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->rowmeans, newsize) );
159  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->rowvariances, newsize) );
160  	
161  	      assert(SCIPgetStage(scip) == SCIP_STAGE_SOLVING);
162  	
163  	      vars = SCIPgetVars(scip);
164  	      nvars = SCIPgetNVars(scip);
165  	
166  	      assert(nvars > 0);
167  	
168  	      /* allocate variable update event processing array storage */
169  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->varfilterposs, nvars) );
170  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->varposs, nvars) );
171  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->updatedvars, nvars) );
172  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->currentubs, nvars) );
173  	      SCIP_CALL( SCIPallocBufferArray(scip, &heurdata->currentlbs, nvars) );
174  	
175  	      heurdata->varpossmemsize = nvars;
176  	      heurdata->nupdatedvars = 0;
177  	
178  	      /* init variable event processing data */
179  	      for( v = 0; v < nvars; ++v )
180  	      {
181  	         assert(SCIPvarIsActive(vars[v]));
182  	         assert(SCIPvarGetProbindex(vars[v]) == v);
183  	
184  	         /* set up variable events to catch bound changes */
185  	         SCIP_CALL( SCIPcatchVarEvent(scip, vars[v], EVENT_DISTRIBUTION, heurdata->eventhdlr, NULL, &(heurdata->varfilterposs[v])) );
186  	         assert(heurdata->varfilterposs[v] >= 0);
187  	
188  	         heurdata->varposs[v] = -1;
189  	         heurdata->updatedvars[v] = NULL;
190  	         heurdata->currentlbs[v] = SCIP_INVALID;
191  	         heurdata->currentubs[v] = SCIP_INVALID;
192  	      }
193  	   }
194  	   else
195  	   {
196  	      SCIP_CALL( SCIPreallocBufferArray(scip, &heurdata->rowinfinitiesdown, newsize) );
197  	      SCIP_CALL( SCIPreallocBufferArray(scip, &heurdata->rowinfinitiesup, newsize) );
198  	      SCIP_CALL( SCIPreallocBufferArray(scip, &heurdata->rowmeans, newsize) );
199  	      SCIP_CALL( SCIPreallocBufferArray(scip, &heurdata->rowvariances, newsize) );
200  	   }
201  	
202  	   /* loop over extended arrays and invalidate data to trigger initialization of this row when necessary */
203  	   for( r = heurdata->memsize; r < newsize; ++r )
204  	   {
205  	      heurdata->rowmeans[r] = SCIP_INVALID;
206  	      heurdata->rowvariances[r] = SCIP_INVALID;
207  	      heurdata->rowinfinitiesdown[r] = 0;
208  	      heurdata->rowinfinitiesup[r] = 0;
209  	   }
210  	
211  	   /* adjust memsize */
212  	   heurdata->memsize = newsize;
213  	
214  	   return SCIP_OKAY;
215  	}
216  	
217  	/** update the variables current lower and upper bound */
218  	static
219  	void heurdataUpdateCurrentBounds(
220  	   SCIP*                 scip,               /**< SCIP data structure */
221  	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
222  	   SCIP_VAR*             var                 /**< the variable to update current bounds */
223  	   )
224  	{
225  	   int varindex;
226  	   SCIP_Real lblocal;
227  	   SCIP_Real ublocal;
228  	
229  	   assert(var != NULL);
230  	
231  	   varindex = SCIPvarGetProbindex(var);
232  	   assert(0 <= varindex && varindex < heurdata->varpossmemsize);
233  	   lblocal = SCIPvarGetLbLocal(var);
234  	   ublocal = SCIPvarGetUbLocal(var);
235  	
236  	   assert(SCIPisFeasLE(scip, lblocal, ublocal));
237  	
238  	   heurdata->currentlbs[varindex] = lblocal;
239  	   heurdata->currentubs[varindex] = ublocal;
240  	}
241  	
242  	/** calculates the initial mean and variance of the row activity normal distribution.
243  	 *
244  	 *  The mean value \f$ \mu \f$ is given by \f$ \mu = \sum_i=1^n c_i * (lb_i +ub_i) / 2 \f$ where
245  	 *  \f$n \f$ is the number of variables, and \f$ c_i, lb_i, ub_i \f$ are the variable coefficient and
246  	 *  bounds, respectively. With the same notation, the variance \f$ \sigma^2 \f$ is given by
247  	 *  \f$ \sigma^2 = \sum_i=1^n c_i^2 * \sigma^2_i \f$, with the variance being
248  	 *  \f$ \sigma^2_i = ((ub_i - lb_i + 1)^2 - 1) / 12 \f$ for integer variables and
249  	 *  \f$ \sigma^2_i = (ub_i - lb_i)^2 / 12 \f$ for continuous variables.
250  	 */
251  	static
252  	void rowCalculateGauss(
253  	   SCIP*                 scip,               /**< SCIP data structure */
254  	   SCIP_HEURDATA*        heurdata,           /**< the heuristic rule data */
255  	   SCIP_ROW*             row,                /**< the row for which the gaussian normal distribution has to be calculated */
256  	   SCIP_Real*            mu,                 /**< pointer to store the mean value of the gaussian normal distribution */
257  	   SCIP_Real*            sigma2,             /**< pointer to store the variance value of the gaussian normal distribution */
258  	   int*                  rowinfinitiesdown,  /**< pointer to store the number of variables with infinite bounds to DECREASE activity */
259  	   int*                  rowinfinitiesup     /**< pointer to store the number of variables with infinite bounds to INCREASE activity */
260  	   )
261  	{
262  	   SCIP_COL** rowcols;
263  	   SCIP_Real* rowvals;
264  	   int nrowvals;
265  	   int c;
266  	
267  	   assert(scip != NULL);
268  	   assert(row != NULL);
269  	   assert(mu != NULL);
270  	   assert(sigma2 != NULL);
271  	   assert(rowinfinitiesup != NULL);
272  	   assert(rowinfinitiesdown != NULL);
273  	
274  	   rowcols = SCIProwGetCols(row);
275  	   rowvals = SCIProwGetVals(row);
276  	   nrowvals = SCIProwGetNNonz(row);
277  	
278  	   assert(nrowvals == 0 || rowcols != NULL);
279  	   assert(nrowvals == 0 || rowvals != NULL);
280  	
281  	   *mu = SCIProwGetConstant(row);
282  	   *sigma2 = 0.0;
283  	   *rowinfinitiesdown = 0;
284  	   *rowinfinitiesup = 0;
285  	
286  	   /* loop over nonzero row coefficients and sum up the variable contributions to mu and sigma2 */
287  	   for( c = 0; c < nrowvals; ++c )
288  	   {
289  	      SCIP_VAR* colvar;
290  	      SCIP_Real colval;
291  	      SCIP_Real colvarlb;
292  	      SCIP_Real colvarub;
293  	      SCIP_Real squarecoeff;
294  	      SCIP_Real varvariance;
295  	      SCIP_Real varmean;
296  	      int varindex;
297  	
298  	      assert(rowcols[c] != NULL);
299  	      colvar = SCIPcolGetVar(rowcols[c]);
300  	      assert(colvar != NULL);
301  	
302  	      colval = rowvals[c];
303  	      colvarlb = SCIPvarGetLbLocal(colvar);
304  	      colvarub = SCIPvarGetUbLocal(colvar);
305  	
306  	      varmean = 0.0;
307  	      varvariance = 0.0;
308  	      varindex = SCIPvarGetProbindex(colvar);
309  	      assert((heurdata->currentlbs[varindex] == SCIP_INVALID) == (heurdata->currentubs[varindex] == SCIP_INVALID)); /*lint !e777 doesn't like comparing floats for equality */
310  	
311  	      /* variable bounds need to be watched from now on */
312  	      if( heurdata->currentlbs[varindex] == SCIP_INVALID ) /*lint !e777 doesn't like comparing floats for equality */
313  	         heurdataUpdateCurrentBounds(scip, heurdata, colvar);
314  	
315  	      assert(!SCIPisInfinity(scip, colvarlb));
316  	      assert(!SCIPisInfinity(scip, -colvarub));
317  	      assert(SCIPisFeasLE(scip, colvarlb, colvarub));
318  	
319  	      /* variables with infinite bounds are skipped for the calculation of the variance; they need to
320  	       * be accounted for by the counters for infinite row activity decrease and increase and they
321  	       * are used to shift the row activity mean in case they have one nonzero, but finite bound */
322  	      if( SCIPisInfinity(scip, -colvarlb) || SCIPisInfinity(scip, colvarub) )
323  	      {
324  	         if( SCIPisInfinity(scip, colvarub) )
325  	         {
326  	         /* an infinite upper bound gives the row an infinite maximum activity or minimum activity, if the coefficient is
327  	          * positive or negative, resp.
328  	          */
329  	            if( colval < 0.0 )
330  	               ++(*rowinfinitiesdown);
331  	            else
332  	               ++(*rowinfinitiesup);
333  	         }
334  	
335  	         /* an infinite lower bound gives the row an infinite maximum activity or minimum activity, if the coefficient is
336  	          * negative or positive, resp.
337  	          */
338  	         if( SCIPisInfinity(scip, -colvarlb) )
339  	         {
340  	            if( colval > 0.0 )
341  	               ++(*rowinfinitiesdown);
342  	            else
343  	               ++(*rowinfinitiesup);
344  	         }
345  	      }
346  	      SCIPvarCalcDistributionParameters(scip, colvarlb, colvarub, SCIPvarGetType(colvar), &varmean, &varvariance);
347  	
348  	      /* actual values are updated; the contribution of the variable to mu is the arithmetic mean of its bounds */
349  	      *mu += colval * varmean;
350  	
351  	      /* the variance contribution of a variable is c^2 * (u - l)^2 / 12.0 for continuous and c^2 * ((u - l + 1)^2 - 1) / 12.0 for integer */
352  	      squarecoeff = SQUARED(colval);
353  	      *sigma2 += squarecoeff * varvariance;
354  	
355  	      assert(!SCIPisFeasNegative(scip, *sigma2));
356  	   }
357  	
358  	   SCIPdebug( SCIPprintRow(scip, row, NULL) );
359  	   SCIPdebugMsg(scip, "  Row %s has a mean value of %g at a sigma2 of %g \n", SCIProwGetName(row), *mu, *sigma2);
360  	}
361  	
362  	/** calculate the branching score of a variable, depending on the chosen score parameter */
363  	static
364  	SCIP_RETCODE calcBranchScore(
365  	   SCIP*                 scip,               /**< current SCIP */
366  	   SCIP_HEURDATA*        heurdata,           /**< branch rule data */
367  	   SCIP_VAR*             var,                /**< candidate variable */
368  	   SCIP_Real             lpsolval,           /**< current fractional LP-relaxation solution value  */
369  	   SCIP_Real*            upscore,            /**< pointer to store the variable score when branching on it in upward direction */
370  	   SCIP_Real*            downscore,          /**< pointer to store the variable score when branching on it in downward direction */
371  	   char                  scoreparam          /**< the score parameter of this heuristic */
372  	   )
373  	{
374  	   SCIP_COL* varcol;
375  	   SCIP_ROW** colrows;
376  	   SCIP_Real* rowvals;
377  	   SCIP_Real varlb;
378  	   SCIP_Real varub;
379  	   SCIP_Real squaredbounddiff; /* current squared difference of variable bounds (ub - lb)^2 */
380  	   SCIP_Real newub;            /* new upper bound if branching downwards */
381  	   SCIP_Real newlb;            /* new lower bound if branching upwards */
382  	   SCIP_Real squaredbounddiffup; /* squared difference after branching upwards (ub - lb')^2 */
383  	   SCIP_Real squaredbounddiffdown; /* squared difference after branching downwards (ub' - lb)^2 */
384  	   SCIP_Real currentmean;      /* current mean value of variable uniform distribution */
385  	   SCIP_Real meanup;           /* mean value of variable uniform distribution after branching up */
386  	   SCIP_Real meandown;         /* mean value of variable uniform distribution after branching down*/
387  	   SCIP_VARTYPE vartype;
388  	   int ncolrows;
389  	   int i;
390  	
391  	   SCIP_Bool onlyactiverows; /* should only rows which are active at the current node be considered? */
392  	
393  	   assert(scip != NULL);
394  	   assert(var != NULL);
395  	   assert(upscore != NULL);
396  	   assert(downscore != NULL);
397  	   assert(!SCIPisIntegral(scip, lpsolval) || SCIPvarIsBinary(var));
398  	   assert(SCIPvarGetStatus(var) == SCIP_VARSTATUS_COLUMN);
399  	
400  	   varcol = SCIPvarGetCol(var);
401  	   assert(varcol != NULL);
402  	
403  	   colrows = SCIPcolGetRows(varcol);
404  	   rowvals = SCIPcolGetVals(varcol);
405  	   ncolrows = SCIPcolGetNNonz(varcol);
406  	   varlb = SCIPvarGetLbLocal(var);
407  	   varub = SCIPvarGetUbLocal(var);
408  	   assert(varub - varlb > 0.5);
409  	   vartype = SCIPvarGetType(var);
410  	
411  	   /* calculate mean and variance of variable uniform distribution before and after branching */
412  	   currentmean = 0.0;
413  	   squaredbounddiff = 0.0;
414  	   SCIPvarCalcDistributionParameters(scip, varlb, varub, vartype, &currentmean, &squaredbounddiff);
415  	
416  	   /* unfixed binary variables may have an integer solution value in the LP solution, eg, at the presence of indicator constraints */
417  	   if( !SCIPvarIsBinary(var) )
418  	   {
419  	      newlb = SCIPfeasCeil(scip, lpsolval);
420  	      newub = SCIPfeasFloor(scip, lpsolval);
421  	   }
422  	   else
423  	   {
424  	      newlb = 1.0;
425  	      newub = 0.0;
426  	   }
427  	
428  	   /* calculate the variable's uniform distribution after branching up and down, respectively. */
429  	   squaredbounddiffup = 0.0;
430  	   meanup = 0.0;
431  	   SCIPvarCalcDistributionParameters(scip, newlb, varub, vartype, &meanup, &squaredbounddiffup);
432  	
433  	   /* calculate the distribution mean and variance for a variable with finite lower bound */
434  	   squaredbounddiffdown = 0.0;
435  	   meandown = 0.0;
436  	   SCIPvarCalcDistributionParameters(scip, varlb, newub, vartype, &meandown, &squaredbounddiffdown);
437  	
438  	   /* initialize the variable's up and down score */
439  	   *upscore = 0.0;
440  	   *downscore = 0.0;
441  	
442  	   onlyactiverows = FALSE;
443  	
444  	   /* loop over the variable rows and calculate the up and down score */
445  	   for( i = 0; i < ncolrows; ++i )
446  	   {
447  	      SCIP_ROW* row;
448  	      SCIP_Real changedrowmean;
449  	      SCIP_Real rowmean;
450  	      SCIP_Real rowvariance;
451  	      SCIP_Real changedrowvariance;
452  	      SCIP_Real currentrowprob;
453  	      SCIP_Real newrowprobup;
454  	      SCIP_Real newrowprobdown;
455  	      SCIP_Real squaredcoeff;
456  	      SCIP_Real rowval;
457  	      int rowinfinitiesdown;
458  	      int rowinfinitiesup;
459  	      int rowpos;
460  	
461  	      row = colrows[i];
462  	      rowval = rowvals[i];
463  	      assert(row != NULL);
464  	
465  	      /* we access the rows by their index */
466  	      rowpos = SCIProwGetIndex(row);
467  	
468  	      /* skip non-active rows if the user parameter was set this way */
469  	      if( onlyactiverows && SCIPisSumPositive(scip, SCIPgetRowLPFeasibility(scip, row)) )
470  	         continue;
471  	
472  	      /* call method to ensure sufficient data capacity */
473  	      SCIP_CALL( heurdataEnsureArraySize(scip, heurdata, rowpos) );
474  	
475  	      /* calculate row activity distribution if this is the first candidate to appear in this row */
476  	      if( heurdata->rowmeans[rowpos] == SCIP_INVALID ) /*lint !e777 doesn't like comparing floats for equality */
477  	      {
478  	         rowCalculateGauss(scip, heurdata, row, &heurdata->rowmeans[rowpos], &heurdata->rowvariances[rowpos],
479  	               &heurdata->rowinfinitiesdown[rowpos], &heurdata->rowinfinitiesup[rowpos]);
480  	      }
481  	
482  	      /* retrieve the row distribution parameters from the branch rule data */
483  	      rowmean = heurdata->rowmeans[rowpos];
484  	      rowvariance = heurdata->rowvariances[rowpos];
485  	      rowinfinitiesdown = heurdata->rowinfinitiesdown[rowpos];
486  	      rowinfinitiesup = heurdata->rowinfinitiesup[rowpos];
487  	      assert(!SCIPisNegative(scip, rowvariance));
488  	
489  	      currentrowprob = SCIProwCalcProbability(scip, row, rowmean, rowvariance,
490  	            rowinfinitiesdown, rowinfinitiesup);
491  	
492  	      /* get variable's current expected contribution to row activity */
493  	      squaredcoeff = SQUARED(rowval);
494  	
495  	      /* first, get the probability change for the row if the variable is branched on upwards. The probability
496  	       * can only be affected if the variable upper bound is finite
497  	       */
498  	      if( !SCIPisInfinity(scip, varub) )
499  	      {
500  	         int rowinftiesdownafterbranch;
501  	         int rowinftiesupafterbranch;
502  	
503  	         /* calculate how branching would affect the row parameters */
504  	         changedrowmean = rowmean + rowval * (meanup - currentmean);
505  	         changedrowvariance = rowvariance + squaredcoeff * (squaredbounddiffup - squaredbounddiff);
506  	         changedrowvariance = MAX(0.0, changedrowvariance);
507  	
508  	         rowinftiesdownafterbranch = rowinfinitiesdown;
509  	         rowinftiesupafterbranch = rowinfinitiesup;
510  	
511  	         /* account for changes of the row's infinite bound contributions */
512  	         if( SCIPisInfinity(scip, -varlb) && rowval < 0.0 )
513  	            rowinftiesupafterbranch--;
514  	         if( SCIPisInfinity(scip, -varlb) && rowval > 0.0 )
515  	            rowinftiesdownafterbranch--;
516  	
517  	         assert(rowinftiesupafterbranch >= 0);
518  	         assert(rowinftiesdownafterbranch >= 0);
519  	         newrowprobup = SCIProwCalcProbability(scip, row, changedrowmean, changedrowvariance, rowinftiesdownafterbranch,
520  	               rowinftiesupafterbranch);
521  	      }
522  	      else
523  	         newrowprobup = currentrowprob;
524  	
525  	      /* do the same for the other branching direction */
526  	      if( !SCIPisInfinity(scip, varlb) )
527  	      {
528  	         int rowinftiesdownafterbranch;
529  	         int rowinftiesupafterbranch;
530  	
531  	         changedrowmean = rowmean + rowval * (meandown - currentmean);
532  	         changedrowvariance = rowvariance + squaredcoeff * (squaredbounddiffdown - squaredbounddiff);
533  	         changedrowvariance = MAX(0.0, changedrowvariance);
534  	
535  	         rowinftiesdownafterbranch = rowinfinitiesdown;
536  	         rowinftiesupafterbranch = rowinfinitiesup;
537  	
538  	         /* account for changes of the row's infinite bound contributions */
539  	         if( SCIPisInfinity(scip, varub) && rowval > 0.0 )
540  	            rowinftiesupafterbranch -= 1;
541  	         if( SCIPisInfinity(scip, varub) && rowval < 0.0 )
542  	            rowinftiesdownafterbranch -= 1;
543  	
544  	         assert(rowinftiesdownafterbranch >= 0);
545  	         assert(rowinftiesupafterbranch >= 0);
546  	         newrowprobdown = SCIProwCalcProbability(scip, row, changedrowmean, changedrowvariance, rowinftiesdownafterbranch,
547  	               rowinftiesupafterbranch);
548  	      }
549  	      else
550  	         newrowprobdown = currentrowprob;
551  	
552  	      /* update the up and down score depending on the chosen scoring parameter */
553  	      SCIP_CALL( SCIPupdateDistributionScore(scip, currentrowprob, newrowprobup, newrowprobdown, upscore, downscore, scoreparam) );
554  	
555  	      SCIPdebugMsg(scip, "  Variable %s changes probability of row %s from %g to %g (branch up) or %g;\n",
556  	         SCIPvarGetName(var), SCIProwGetName(row), currentrowprob, newrowprobup, newrowprobdown);
557  	      SCIPdebugMsg(scip, "  -->  new variable score: %g (for branching up), %g (for branching down)\n",
558  	         *upscore, *downscore);
559  	   }
560  	
561  	   return SCIP_OKAY;
562  	}
563  	
564  	/** free heuristic data */
565  	static
566  	SCIP_RETCODE heurdataFreeArrays(
567  	   SCIP*                 scip,               /**< SCIP data structure */
568  	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
569  	   )
570  	{
571  	   assert(heurdata->memsize == 0 || heurdata->rowmeans != NULL);
572  	   assert(heurdata->memsize >= 0);
573  	
574  	   if( heurdata->varpossmemsize > 0 )
575  	   {
576  	      SCIP_VAR** vars;
577  	      int v;
578  	
579  	      assert(heurdata->varpossmemsize == SCIPgetNVars(scip));
580  	
581  	      vars = SCIPgetVars(scip);
582  	      for( v = heurdata->varpossmemsize - 1; v >= 0; --v )
583  	      {
584  	         SCIP_VAR* var;
585  	
586  	         var = vars[v];
587  	
588  	         assert(var != NULL);
589  	         assert(v == SCIPvarGetProbindex(var));
590  	         SCIP_CALL( SCIPdropVarEvent(scip, var, EVENT_DISTRIBUTION, heurdata->eventhdlr, NULL, heurdata->varfilterposs[v]) );
591  	      }
592  	      SCIPfreeBufferArray(scip, &heurdata->currentlbs);
593  	      SCIPfreeBufferArray(scip, &heurdata->currentubs);
594  	      SCIPfreeBufferArray(scip, &heurdata->updatedvars);
595  	      SCIPfreeBufferArray(scip, &heurdata->varposs);
596  	      SCIPfreeBufferArray(scip, &heurdata->varfilterposs);
597  	   }
598  	
599  	   if( heurdata->memsize > 0 )
600  	   {
601  	      SCIPfreeBufferArray(scip, &heurdata->rowvariances);
602  	      SCIPfreeBufferArray(scip, &heurdata->rowmeans);
603  	      SCIPfreeBufferArray(scip, &heurdata->rowinfinitiesup);
604  	      SCIPfreeBufferArray(scip, &heurdata->rowinfinitiesdown);
605  	
606  	      heurdata->memsize = 0;
607  	   }
608  	
609  	   heurdata->varpossmemsize = 0;
610  	   heurdata->nupdatedvars = 0;
611  	
612  	   return SCIP_OKAY;
613  	}
614  	
615  	/** add variable to the bound change event queue; skipped if variable is already in there, or if variable has
616  	 *  no row currently watched
617  	 */
618  	static
619  	void heurdataAddBoundChangeVar(
620  	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
621  	   SCIP_VAR*             var                 /**< the variable whose bound changes need to be processed */
622  	   )
623  	{
624  	   int varindex;
625  	   int varpos;
626  	
627  	   assert(var != NULL);
628  	
629  	   varindex = SCIPvarGetProbindex(var);
630  	   assert(-1 <= varindex && varindex < heurdata->varpossmemsize);
631  	
632  	   /* if variable is not active, it should not be watched */
633  	   if( varindex == -1 )
634  	      return;
635  	   varpos = heurdata->varposs[varindex];
636  	   assert(varpos < heurdata->nupdatedvars);
637  	
638  	   /* nothing to do if variable is already in the queue */
639  	   if( varpos >= 0 )
640  	   {
641  	      assert(heurdata->updatedvars[varpos] == var);
642  	
643  	      return;
644  	   }
645  	
646  	   /* if none of the variables rows was calculated yet, variable needs not to be watched */
647  	   assert((heurdata->currentlbs[varindex] == SCIP_INVALID) == (heurdata->currentubs[varindex] == SCIP_INVALID)); /*lint !e777 doesn't like comparing floats for equality */
648  	
649  	   /* we don't need to enqueue the variable if it hasn't been watched so far */
650  	   if( heurdata->currentlbs[varindex] == SCIP_INVALID ) /*lint !e777 see above */
651  	      return;
652  	
653  	   /* add the variable to the heuristic data of variables to process updates for */
654  	   assert(heurdata->varpossmemsize > heurdata->nupdatedvars);
655  	   varpos = heurdata->nupdatedvars;
656  	   heurdata->updatedvars[varpos] = var;
657  	   heurdata->varposs[varindex] = varpos;
658  	   ++heurdata->nupdatedvars;
659  	}
660  	
661  	/** returns the next unprocessed variable (last in, first out) with pending bound changes, or NULL */
662  	static
663  	SCIP_VAR* heurdataPopBoundChangeVar(
664  	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
665  	   )
666  	{
667  	   SCIP_VAR* var;
668  	   int varpos;
669  	   int varindex;
670  	
671  	   assert(heurdata->nupdatedvars >= 0);
672  	
673  	   /* return if no variable is currently pending */
674  	   if( heurdata->nupdatedvars == 0 )
675  	      return NULL;
676  	
677  	   varpos = heurdata->nupdatedvars - 1;
678  	   var = heurdata->updatedvars[varpos];
679  	   assert(var != NULL);
680  	   varindex = SCIPvarGetProbindex(var);
681  	   assert(0 <= varindex && varindex < heurdata->varpossmemsize);
682  	   assert(varpos == heurdata->varposs[varindex]);
683  	
684  	   heurdata->varposs[varindex] = -1;
685  	   heurdata->nupdatedvars--;
686  	
687  	   return var;
688  	}
689  	
690  	/** process a variable from the queue of changed variables */
691  	static
692  	SCIP_RETCODE varProcessBoundChanges(
693  	   SCIP*                 scip,               /**< SCIP data structure */
694  	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
695  	   SCIP_VAR*             var                 /**< the variable whose bound changes need to be processed */
696  	   )
697  	{
698  	   SCIP_ROW** colrows;
699  	   SCIP_COL* varcol;
700  	   SCIP_Real* colvals;
701  	   SCIP_Real oldmean;
702  	   SCIP_Real newmean;
703  	   SCIP_Real oldvariance;
704  	   SCIP_Real newvariance;
705  	   SCIP_Real oldlb;
706  	   SCIP_Real newlb;
707  	   SCIP_Real oldub;
708  	   SCIP_Real newub;
709  	   SCIP_VARTYPE vartype;
710  	   int ncolrows;
711  	   int r;
712  	   int varindex;
713  	
714  	   /* ensure that this is a probing bound change */
715  	   assert(SCIPinProbing(scip));
716  	
717  	   assert(var != NULL);
718  	   varcol = SCIPvarGetCol(var);
719  	   assert(varcol != NULL);
720  	   colrows = SCIPcolGetRows(varcol);
721  	   colvals = SCIPcolGetVals(varcol);
722  	   ncolrows = SCIPcolGetNNonz(varcol);
723  	
724  	   varindex = SCIPvarGetProbindex(var);
725  	
726  	   oldlb = heurdata->currentlbs[varindex];
727  	   oldub = heurdata->currentubs[varindex];
728  	
729  	   /* skip update if the variable has never been subject of previously calculated row activities */
730  	   assert((oldlb == SCIP_INVALID) == (oldub == SCIP_INVALID)); /*lint !e777 doesn't like comparing floats for equality */
731  	   if( oldlb == SCIP_INVALID ) /*lint !e777 */
732  	      return SCIP_OKAY;
733  	
734  	   newlb = SCIPvarGetLbLocal(var);
735  	   newub = SCIPvarGetUbLocal(var);
736  	
737  	   /* skip update if the bound change events have cancelled out */
738  	   if( SCIPisFeasEQ(scip, oldlb, newlb) && SCIPisFeasEQ(scip, oldub, newub) )
739  	      return SCIP_OKAY;
740  	
741  	   /* calculate old and new variable distribution mean and variance */
742  	   oldvariance = 0.0;
743  	   newvariance = 0.0;
744  	   oldmean = 0.0;
745  	   newmean = 0.0;
746  	   vartype = SCIPvarGetType(var);
747  	   SCIPvarCalcDistributionParameters(scip, oldlb, oldub, vartype, &oldmean, &oldvariance);
748  	   SCIPvarCalcDistributionParameters(scip, newlb, newub, vartype, &newmean, &newvariance);
749  	
750  	   /* loop over all rows of this variable and update activity distribution */
751  	   for( r = 0; r < ncolrows; ++r )
752  	   {
753  	      int rowpos;
754  	
755  	      assert(colrows[r] != NULL);
756  	      rowpos = SCIProwGetIndex(colrows[r]);
757  	      assert(rowpos >= 0);
758  	
759  	      SCIP_CALL( heurdataEnsureArraySize(scip, heurdata, rowpos) );
760  	
761  	      /* only consider rows for which activity distribution was already calculated */
762  	      if( heurdata->rowmeans[rowpos] != SCIP_INVALID ) /*lint !e777 doesn't like comparing floats for equality */
763  	      {
764  	         SCIP_Real coeff;
765  	         SCIP_Real coeffsquared;
766  	         assert(heurdata->rowvariances[rowpos] != SCIP_INVALID && SCIPisFeasGE(scip, heurdata->rowvariances[rowpos], 0.0)); /*lint !e777 */
767  	
768  	         coeff = colvals[r];
769  	         coeffsquared = SQUARED(coeff);
770  	
771  	         /* update variable contribution to row activity distribution */
772  	         heurdata->rowmeans[rowpos] += coeff * (newmean - oldmean);
773  	         heurdata->rowvariances[rowpos] += coeffsquared * (newvariance - oldvariance);
774  	         heurdata->rowvariances[rowpos] = MAX(0.0, heurdata->rowvariances[rowpos]);
775  	
776  	         /* account for changes of the infinite contributions to row activities */
777  	         if( coeff > 0.0 )
778  	         {
779  	            /* if the coefficient is positive, upper bounds affect activity up */
780  	            if( SCIPisInfinity(scip, newub) && !SCIPisInfinity(scip, oldub) )
781  	               ++heurdata->rowinfinitiesup[rowpos];
782  	            else if( !SCIPisInfinity(scip, newub) && SCIPisInfinity(scip, oldub) )
783  	               --heurdata->rowinfinitiesup[rowpos];
784  	
785  	            if( SCIPisInfinity(scip, newlb) && !SCIPisInfinity(scip, oldlb) )
786  	               ++heurdata->rowinfinitiesdown[rowpos];
787  	            else if( !SCIPisInfinity(scip, newlb) && SCIPisInfinity(scip, oldlb) )
788  	               --heurdata->rowinfinitiesdown[rowpos];
789  	         }
790  	         else if( coeff < 0.0 )
791  	         {
792  	            if( SCIPisInfinity(scip, newub) && !SCIPisInfinity(scip, oldub) )
793  	               ++heurdata->rowinfinitiesdown[rowpos];
794  	            else if( !SCIPisInfinity(scip, newub) && SCIPisInfinity(scip, oldub) )
795  	               --heurdata->rowinfinitiesdown[rowpos];
796  	
797  	            if( SCIPisInfinity(scip, newlb) && !SCIPisInfinity(scip, oldlb) )
798  	               ++heurdata->rowinfinitiesup[rowpos];
799  	            else if( !SCIPisInfinity(scip, newlb) && SCIPisInfinity(scip, oldlb) )
800  	               --heurdata->rowinfinitiesup[rowpos];
801  	         }
802  	         assert(heurdata->rowinfinitiesdown[rowpos] >= 0);
803  	         assert(heurdata->rowinfinitiesup[rowpos] >= 0);
804  	      }
805  	   }
806  	
807  	   /* store the new local bounds in the data */
808  	   heurdataUpdateCurrentBounds(scip, heurdata, var);
809  	
810  	   return SCIP_OKAY;
811  	}
812  	
813  	/** destructor of event handler to free user data (called when SCIP is exiting) */
814  	static
815  	SCIP_DECL_EVENTFREE(eventFreeDistributiondiving)
816  	{
817  	   SCIP_EVENTHDLRDATA* eventhdlrdata;
818  	
819  	   eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
820  	   assert(eventhdlrdata != NULL);
821  	
822  	   SCIPfreeBlockMemory(scip, &eventhdlrdata);
823  	   SCIPeventhdlrSetData(eventhdlr, NULL);
824  	
825  	   return SCIP_OKAY;
826  	}
827  	
828  	/*
829  	 * Callback methods
830  	 */
831  	
832  	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
833  	static
834  	SCIP_DECL_HEURCOPY(heurCopyDistributiondiving)
835  	{  /*lint --e{715}*/
836  	   assert(scip != NULL);
837  	   assert(heur != NULL);
838  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
839  	
840  	   /* call inclusion method of primal heuristic */
841  	   SCIP_CALL( SCIPincludeHeurDistributiondiving(scip) );
842  	
843  	   return SCIP_OKAY;
844  	}
845  	
846  	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
847  	static
848  	SCIP_DECL_HEURFREE(heurFreeDistributiondiving) /*lint --e{715}*/
849  	{  /*lint --e{715}*/
850  	   SCIP_HEURDATA* heurdata;
851  	
852  	   assert(heur != NULL);
853  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
854  	   assert(scip != NULL);
855  	
856  	   /* free heuristic data */
857  	   heurdata = SCIPheurGetData(heur);
858  	   assert(heurdata != NULL);
859  	   SCIPfreeBlockMemory(scip, &heurdata);
860  	   SCIPheurSetData(heur, NULL);
861  	
862  	   return SCIP_OKAY;
863  	}
864  	
865  	
866  	/** initialization method of primal heuristic (called after problem was transformed) */
867  	static
868  	SCIP_DECL_HEURINIT(heurInitDistributiondiving) /*lint --e{715}*/
869  	{  /*lint --e{715}*/
870  	   SCIP_HEURDATA* heurdata;
871  	
872  	   assert(heur != NULL);
873  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
874  	
875  	   /* get heuristic data */
876  	   heurdata = SCIPheurGetData(heur);
877  	   assert(heurdata != NULL);
878  	
879  	   /* create working solution */
880  	   SCIP_CALL( SCIPcreateSol(scip, &heurdata->sol, heur) );
881  	
882  	   return SCIP_OKAY;
883  	}
884  	
885  	
886  	/** deinitialization method of primal heuristic (called before transformed problem is freed) */
887  	static
888  	SCIP_DECL_HEUREXIT(heurExitDistributiondiving) /*lint --e{715}*/
889  	{  /*lint --e{715}*/
890  	   SCIP_HEURDATA* heurdata;
891  	
892  	   assert(heur != NULL);
893  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
894  	
895  	   /* get heuristic data */
896  	   heurdata = SCIPheurGetData(heur);
897  	   assert(heurdata != NULL);
898  	
899  	   /* free working solution */
900  	   SCIP_CALL( SCIPfreeSol(scip, &heurdata->sol) );
901  	
902  	   return SCIP_OKAY;
903  	}
904  	
905  	/** scoring callback for distribution diving. best candidate maximizes the distribution score */
906  	static
907  	SCIP_DECL_DIVESETGETSCORE(divesetGetScoreDistributiondiving)
908  	{  /*lint --e{715}*/
909  	   SCIP_HEURDATA* heurdata;
910  	   SCIP_Real upscore;
911  	   SCIP_Real downscore;
912  	   int varindex;
913  	
914  	   heurdata = SCIPheurGetData(SCIPdivesetGetHeur(diveset));
915  	   assert(heurdata != NULL);
916  	
917  	   /* process pending bound change events */
918  	   while( heurdata->nupdatedvars > 0 )
919  	   {
920  	      SCIP_VAR* nextvar;
921  	
922  	      /* pop the next variable from the queue and process its bound changes */
923  	      nextvar = heurdataPopBoundChangeVar(heurdata);
924  	      assert(nextvar != NULL);
925  	      SCIP_CALL( varProcessBoundChanges(scip, heurdata, nextvar) );
926  	   }
927  	
928  	   assert(cand != NULL);
929  	
930  	   varindex = SCIPvarGetProbindex(cand);
931  	
932  	   /* terminate with a penalty for inactive variables, which the plugin can currently not score
933  	    * this should never happen with default settings where only LP branching candidates are iterated, but might occur
934  	    * if other constraint handlers try to score an inactive variable that was (multi-)aggregated or negated
935  	    */
936  	   if( varindex == - 1 )
937  	   {
938  	      *score = -1.0;
939  	      *roundup = FALSE;
940  	
941  	      return SCIP_OKAY;
942  	   }
943  	
944  	   /* in debug mode, ensure that all bound process events which occurred in the mean time have been captured
945  	    * by the heuristic event system
946  	    */
947  	   assert(SCIPisFeasLE(scip, SCIPvarGetLbLocal(cand), SCIPvarGetUbLocal(cand)));
948  	   assert(0 <= varindex && varindex < heurdata->varpossmemsize);
949  	
950  	   assert((heurdata->currentlbs[varindex] == SCIP_INVALID) == (heurdata->currentubs[varindex] == SCIP_INVALID));/*lint !e777 doesn't like comparing floats for equality */
951  	   assert((heurdata->currentlbs[varindex] == SCIP_INVALID) || SCIPisFeasEQ(scip, SCIPvarGetLbLocal(cand), heurdata->currentlbs[varindex])); /*lint !e777 */
952  	   assert((heurdata->currentubs[varindex] == SCIP_INVALID) || SCIPisFeasEQ(scip, SCIPvarGetUbLocal(cand), heurdata->currentubs[varindex])); /*lint !e777 */
953  	
954  	   /* if the heuristic has not captured the variable bounds yet, this can be done now */
955  	   if( heurdata->currentlbs[varindex] == SCIP_INVALID ) /*lint !e777 */
956  	      heurdataUpdateCurrentBounds(scip, heurdata, cand);
957  	
958  	   upscore = 0.0;
959  	   downscore = 0.0;
960  	
961  	   /* loop over candidate rows and determine the candidate up- and down- branching score w.r.t. the score parameter */
962  	   SCIP_CALL( calcBranchScore(scip, heurdata, cand, candsol, &upscore, &downscore, heurdata->score) );
963  	
964  	   /* score is simply the maximum of the two individual scores */
965  	   *roundup = (upscore > downscore);
966  	   *score = MAX(upscore, downscore);
967  	
968  	   return SCIP_OKAY;
969  	}
970  	
971  	
972  	/** execution method of primal heuristic */
973  	static
974  	SCIP_DECL_HEUREXEC(heurExecDistributiondiving)
975  	{  /*lint --e{715}*/
976  	   SCIP_HEURDATA* heurdata;
977  	   SCIP_DIVESET* diveset;
978  	   int nlprows;
979  	
980  	   assert(heur != NULL);
981  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
982  	   assert(scip != NULL);
983  	   assert(result != NULL);
984  	   assert(SCIPhasCurrentNodeLP(scip));
985  	
986  	   *result = SCIP_DIDNOTRUN;
987  	
988  	   /* get heuristic's data */
989  	   heurdata = SCIPheurGetData(heur);
990  	   assert(heurdata != NULL);
991  	   nlprows = SCIPgetNLPRows(scip);
992  	   if( nlprows == 0 )
993  	      return SCIP_OKAY;
994  	
995  	   /* terminate if there are no integer variables (note that, e.g., SOS1 variables may be present) */
996  	   if( SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) == 0 )
997  	      return SCIP_OKAY;
998  	
999  	   /* select and store the scoring parameter for this call of the heuristic */
1000 	   if( heurdata->scoreparam == 'r' )
1001 	      heurdata->score = SCOREPARAM_VALUES[SCIPheurGetNCalls(heur) % SCOREPARAM_VALUESLEN];
1002 	   else
1003 	      heurdata->score = heurdata->scoreparam;
1004 	
1005 	   SCIP_CALL( heurdataEnsureArraySize(scip, heurdata, nlprows) );
1006 	   assert(SCIPheurGetNDivesets(heur) > 0);
1007 	   assert(SCIPheurGetDivesets(heur) != NULL);
1008 	   diveset = SCIPheurGetDivesets(heur)[0];
1009 	   assert(diveset != NULL);
1010 	
1011 	   SCIP_CALL( SCIPperformGenericDivingAlgorithm(scip, diveset, heurdata->sol, heur, result, nodeinfeasible, -1L, -1, -1.0, SCIP_DIVECONTEXT_SINGLE) );
1012 	
1013 	   SCIP_CALL( heurdataFreeArrays(scip, heurdata) );
1014 	
1015 	   return SCIP_OKAY;
1016 	}
1017 	
1018 	/** event execution method of distribution branching which handles bound change events of variables */
1019 	static
1020 	SCIP_DECL_EVENTEXEC(eventExecDistribution)
1021 	{
1022 	   SCIP_HEURDATA* heurdata;
1023 	   SCIP_EVENTHDLRDATA* eventhdlrdata;
1024 	   SCIP_VAR* var;
1025 	
1026 	   assert(eventhdlr != NULL);
1027 	   eventhdlrdata = SCIPeventhdlrGetData(eventhdlr);
1028 	   assert(eventhdlrdata != NULL);
1029 	
1030 	   heurdata = eventhdlrdata->heurdata;
1031 	   var = SCIPeventGetVar(event);
1032 	
1033 	   /* add the variable to the queue of unprocessed variables; method itself ensures that every variable is added at most once */
1034 	   heurdataAddBoundChangeVar(heurdata, var);
1035 	
1036 	   return SCIP_OKAY;
1037 	}
1038 	
1039 	/*
1040 	 * heuristic specific interface methods
1041 	 */
1042 	
1043 	#define divesetAvailableDistributiondiving NULL
1044 	
1045 	/** creates the distributiondiving heuristic and includes it in SCIP */
1046 	SCIP_RETCODE SCIPincludeHeurDistributiondiving(
1047 	   SCIP*                 scip                /**< SCIP data structure */
1048 	   )
1049 	{
1050 	   SCIP_HEURDATA* heurdata;
1051 	   SCIP_HEUR* heur;
1052 	   SCIP_EVENTHDLRDATA* eventhdlrdata;
1053 	
1054 	   /* create distributiondiving data */
1055 	   heurdata = NULL;
1056 	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
1057 	
1058 	   heurdata->memsize = 0;
1059 	   heurdata->rowmeans = NULL;
1060 	   heurdata->rowvariances = NULL;
1061 	   heurdata->rowinfinitiesdown = NULL;
1062 	   heurdata->rowinfinitiesup = NULL;
1063 	   heurdata->varfilterposs = NULL;
1064 	   heurdata->currentlbs = NULL;
1065 	   heurdata->currentubs = NULL;
1066 	
1067 	   /* create event handler first to finish heuristic data */
1068 	   eventhdlrdata = NULL;
1069 	   SCIP_CALL( SCIPallocBlockMemory(scip, &eventhdlrdata) );
1070 	   eventhdlrdata->heurdata = heurdata;
1071 	
1072 	   heurdata->eventhdlr = NULL;
1073 	   SCIP_CALL( SCIPincludeEventhdlrBasic(scip, &heurdata->eventhdlr, EVENTHDLR_NAME,
1074 	         "event handler for dynamic acitivity distribution updating",
1075 	         eventExecDistribution, eventhdlrdata) );
1076 	   assert( heurdata->eventhdlr != NULL);
1077 	   SCIP_CALL( SCIPsetEventhdlrFree(scip, heurdata->eventhdlr, eventFreeDistributiondiving) );
1078 	
1079 	   /* include primal heuristic */
1080 	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
1081 	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
1082 	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecDistributiondiving, heurdata) );
1083 	
1084 	   assert(heur != NULL);
1085 	
1086 	   /* set non-NULL pointers to callback methods */
1087 	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyDistributiondiving) );
1088 	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeDistributiondiving) );
1089 	   SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitDistributiondiving) );
1090 	   SCIP_CALL( SCIPsetHeurExit(scip, heur, heurExitDistributiondiving) );
1091 	
1092 	   /* add diveset with the defined scoring function */
1093 	   SCIP_CALL( SCIPcreateDiveset(scip, NULL, heur, HEUR_NAME, DEFAULT_MINRELDEPTH,
1094 	         DEFAULT_MAXRELDEPTH, DEFAULT_MAXLPITERQUOT, DEFAULT_MAXDIVEUBQUOT,
1095 	         DEFAULT_MAXDIVEAVGQUOT, DEFAULT_MAXDIVEUBQUOTNOSOL,
1096 	         DEFAULT_MAXDIVEAVGQUOTNOSOL, DEFAULT_LPRESOLVEDOMCHGQUOT, DEFAULT_LPSOLVEFREQ,
1097 	         DEFAULT_MAXLPITEROFS, DEFAULT_RANDSEED, DEFAULT_BACKTRACK, DEFAULT_ONLYLPBRANCHCANDS, DIVESET_ISPUBLIC, DIVESET_DIVETYPES,
1098 	         divesetGetScoreDistributiondiving, divesetAvailableDistributiondiving) );
1099 	
1100 	   SCIP_CALL( SCIPaddCharParam(scip, "heuristics/" HEUR_NAME "/scoreparam",
1101 	         "the score;largest 'd'ifference, 'l'owest cumulative probability,'h'ighest c.p., 'v'otes lowest c.p., votes highest c.p.('w'), 'r'evolving",
1102 	         &heurdata->scoreparam, TRUE, DEFAULT_SCOREPARAM, "lvdhwr", NULL, NULL) );
1103 	
1104 	   return SCIP_OKAY;
1105 	}
1106