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_octane.c
26   	 * @ingroup DEFPLUGINS_HEUR
27   	 * @brief  octane primal heuristic based on Balas, Ceria, Dawande, Margot, and Pataki
28   	 * @author Timo Berthold
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_octane.h"
35   	#include "scip/pub_heur.h"
36   	#include "scip/pub_lp.h"
37   	#include "scip/pub_message.h"
38   	#include "scip/pub_misc_sort.h"
39   	#include "scip/pub_var.h"
40   	#include "scip/scip_branch.h"
41   	#include "scip/scip_general.h"
42   	#include "scip/scip_heur.h"
43   	#include "scip/scip_lp.h"
44   	#include "scip/scip_mem.h"
45   	#include "scip/scip_message.h"
46   	#include "scip/scip_numerics.h"
47   	#include "scip/scip_param.h"
48   	#include "scip/scip_prob.h"
49   	#include "scip/scip_sol.h"
50   	#include "scip/scip_solvingstats.h"
51   	#include <string.h>
52   	
53   	#define HEUR_NAME             "octane"
54   	#define HEUR_DESC             "octane primal heuristic for pure {0;1}-problems based on Balas et al."
55   	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_ROUNDING
56   	#define HEUR_PRIORITY         -1008000
57   	#define HEUR_FREQ             -1
58   	#define HEUR_FREQOFS          0
59   	#define HEUR_MAXDEPTH         -1
60   	#define HEUR_TIMING           SCIP_HEURTIMING_AFTERLPNODE
61   	#define HEUR_USESSUBSCIP      FALSE  /**< does the heuristic use a secondary SCIP instance? */
62   	
63   	#define DEFAULT_FMAX          100       /**< {0,1}-points to be checked */
64   	#define DEFAULT_FFIRST        10        /**< {0,1}-points to be generated at first */
65   	#define DEFAULT_USEFRACSPACE  TRUE      /**< use heuristic for the space of fractional variables or for whole space? */
66   	
67   	/** primal heuristic data */
68   	struct SCIP_HeurData
69   	{
70   	   SCIP_SOL*             sol;                /**< working solution */
71   	   int                   f_max;              /**< {0,1}-points to be checked */
72   	   int                   f_first;            /**< {0,1}-points to be generated at first in order to check whether restart is necessary */
73   	   int                   lastrule;           /**< last ray selection rule that was performed */
74   	   SCIP_Bool             usefracspace;       /**< use heuristic for the space of fractional variables or for the whole space? */
75   	   SCIP_Bool             useobjray;          /**< should the inner normal of the objective be used as one ray direction? */
76   	   SCIP_Bool             useavgray;          /**< should the average ray of the basic cone be used as one ray direction? */
77   	   SCIP_Bool             usediffray;         /**< should difference between root sol and current LP sol be used as one ray direction? */
78   	   SCIP_Bool             useavgwgtray;       /**< should the weighted average ray of the basic cone be used as one ray direction? */
79   	   SCIP_Bool             useavgnbray;        /**< should the average ray of the nonbasic cone be used as one ray direction? */
80   	   int                   nsuccess;           /**< number of runs that produced at least one feasible solution */
81   	};
82   	
83   	/*
84   	 * Local methods
85   	 */
86   	
87   	
88   	/** tries to insert the facet obtained from facet i flipped in component j into the list of the fmax nearest facets */
89   	static
90   	void tryToInsert(
91   	   SCIP*                 scip,               /**< SCIP data structure                        */
92   	   SCIP_Bool**           facets,             /**< facets got so far                          */
93   	   SCIP_Real*            lambda,             /**< distances of the facets                    */
94   	   int                   i,                  /**< current facet                              */
95   	   int                   j,                  /**< component to flip                          */
96   	   int                   f_max,              /**< maximal number of facets to create         */
97   	   int                   nsubspacevars,      /**< dimension of the fractional space          */
98   	   SCIP_Real             lam,                /**< distance of the current facet              */
99   	   int*                  nfacets             /**< number of facets                           */
100  	   )
101  	{
102  	   SCIP_Bool* lastfacet;
103  	   int k;
104  	
105  	   assert(scip != NULL);
106  	   assert(facets != NULL);
107  	   assert(lambda != NULL);
108  	   assert(nfacets != NULL);
109  	
110  	   if( SCIPisFeasLE(scip, lam, 0.0) || SCIPisFeasGE(scip, lam, lambda[f_max-1]) )
111  	      return;
112  	
113  	   lastfacet = facets[f_max];
114  	
115  	   /* shifting lam through lambda, lambda keeps increasingly sorted */
116  	   for( k = f_max; k > 0 && SCIPisFeasGT(scip, lambda[k-1], lam); --k )
117  	   {
118  	      lambda[k] = lambda[k-1];
119  	      facets[k] = facets[k-1];
120  	   }
121  	   assert(i < k && k < f_max );
122  	
123  	   /* inserting new facet into list, new facet is facet at position i flipped in coordinate j, new distance lam */
124  	   facets[k] = lastfacet;
125  	   lambda[k] = lam;
126  	
127  	   /*lint --e{866}*/
128  	   BMScopyMemoryArray(facets[k], facets[i], nsubspacevars);
129  	   facets[k][j] = !facets[k][j];
130  	   (*nfacets)++;
131  	}
132  	
133  	/** constructs a solution from a given facet paying attention to the transformations made at the beginning of OCTANE */
134  	static
135  	SCIP_RETCODE getSolFromFacet(
136  	   SCIP*                 scip,               /**< SCIP data structure                   */
137  	   SCIP_Bool*            facet,              /**< current facet                         */
138  	   SCIP_SOL*             sol,                /**< solution to create                    */
139  	   SCIP_Bool*            sign,               /**< marker for retransformation            */
140  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
141  	   int                   nsubspacevars       /**< dimension of fractional space         */
142  	   )
143  	{
144  	   int v;
145  	
146  	   assert(scip != NULL);
147  	   assert(facet != NULL);
148  	   assert(sol != NULL);
149  	   assert(sign != NULL);
150  	   assert(subspacevars != NULL);
151  	
152  	   SCIP_CALL( SCIPlinkLPSol(scip, sol) );
153  	   for( v = nsubspacevars - 1; v >= 0; --v )
154  	   {
155  	      /* after permutation, a variable should be set to 1, iff there was no reflection in this coordinate and the hit
156  	       * facet has coordinate + or there was a reflection and the facet has coordinate - */
157  	      if( facet[v] == sign[v] )
158  	      {
159  	         SCIP_CALL( SCIPsetSolVal(scip, sol, subspacevars[v], 1.0) );
160  	      }
161  	      else
162  	      {
163  	         SCIP_CALL( SCIPsetSolVal(scip, sol, subspacevars[v], 0.0) );
164  	      }
165  	   }
166  	
167  	   return SCIP_OKAY;
168  	}
169  	
170  	/** generates the direction of the shooting ray as the inner normal of the objective function */
171  	static
172  	SCIP_RETCODE generateObjectiveRay(
173  	   SCIP*                 scip,               /**< SCIP data structure                   */
174  	   SCIP_Real*            raydirection,       /**< shooting ray                          */
175  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
176  	   int                   nsubspacevars       /**< dimension of fractional space         */
177  	   )
178  	{
179  	   int v;
180  	
181  	   assert(scip != NULL);
182  	   assert(raydirection != NULL);
183  	   assert(subspacevars != NULL);
184  	
185  	   for( v = nsubspacevars - 1; v >= 0; --v )
186  	      raydirection[v] = SCIPvarGetObj(subspacevars[v]);
187  	   return SCIP_OKAY;
188  	}
189  	
190  	/** generates the direction of the shooting ray as the difference between the root and the current LP solution */
191  	static
192  	SCIP_RETCODE generateDifferenceRay(
193  	   SCIP*                 scip,               /**< SCIP data structure                   */
194  	   SCIP_Real*            raydirection,       /**< shooting ray                          */
195  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
196  	   int                   nsubspacevars       /**< dimension of fractional space         */
197  	   )
198  	{
199  	   int v;
200  	
201  	   assert(scip != NULL);
202  	   assert(raydirection != NULL);
203  	   assert(subspacevars != NULL);
204  	
205  	   for( v = nsubspacevars - 1; v >= 0; --v )
206  	      raydirection[v] = SCIPvarGetLPSol(subspacevars[v]) - SCIPvarGetRootSol(subspacevars[v]);
207  	
208  	   return SCIP_OKAY;
209  	}
210  	
211  	
212  	/** generates the direction of the shooting ray as the average of the extreme rays of the basic cone */
213  	static
214  	SCIP_RETCODE generateAverageRay(
215  	   SCIP*                 scip,               /**< SCIP data structure                   */
216  	   SCIP_Real*            raydirection,       /**< shooting ray                          */
217  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
218  	   int                   nsubspacevars,      /**< dimension of fractional space         */
219  	   SCIP_Bool             weighted            /**< should the rays be weighted?          */
220  	   )
221  	{
222  	   SCIP_ROW** rows;
223  	   SCIP_Real** tableaurows;
224  	   SCIP_Real* rownorm;
225  	   SCIP_Real rowweight;
226  	   int** tableaurowinds;                     /* indices of non-zero entries */
227  	   int* ntableaurowinds;                     /* number of non-zero entries */
228  	   SCIP_Bool* usedrowids = NULL;             /* visited row indices */
229  	   int* rowids;                              /* row indices */
230  	   int nrowids = 0;                          /* number of row indices */
231  	   int tableaurowind;
232  	   int nrows;
233  	   int i;
234  	   int j;
235  	   int sparse = -1; /* used to check that either all information is sparse or not sparse */
236  	
237  	   assert(scip != NULL);
238  	   assert(raydirection != NULL);
239  	   assert(subspacevars != NULL);
240  	   assert(nsubspacevars > 0);
241  	
242  	   /* get data */
243  	   SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
244  	   assert(nrows > 0);
245  	   assert(rows != NULL);
246  	
247  	   /* allocate memory */
248  	   SCIP_CALL( SCIPallocBufferArray(scip, &tableaurows, nsubspacevars) );
249  	   SCIP_CALL( SCIPallocBufferArray(scip, &tableaurowinds, nsubspacevars) );
250  	   SCIP_CALL( SCIPallocBufferArray(scip, &ntableaurowinds, nsubspacevars) );
251  	   for( j = nsubspacevars - 1; j >= 0; --j )
252  	   {
253  	      /*lint --e{866}*/
254  	      SCIP_CALL( SCIPallocBufferArray(scip, &tableaurows[j], nrows) );
255  	      SCIP_CALL( SCIPallocBufferArray(scip, &tableaurowinds[j], nrows) );
256  	   }
257  	
258  	   SCIP_CALL( SCIPallocBufferArray(scip, &rownorm, nrows) );
259  	   BMSclearMemoryArray(rownorm, nrows);
260  	
261  	   /* clear ray */
262  	   BMSclearMemoryArray(raydirection, nsubspacevars);
263  	
264  	   /* get the relevant columns of the simplex tableau */
265  	   for( j = nsubspacevars - 1; j >= 0; --j )
266  	   {
267  	      assert(SCIPcolGetLPPos(SCIPvarGetCol(subspacevars[j])) >= 0);
268  	      SCIP_CALL( SCIPgetLPBInvACol(scip, SCIPcolGetLPPos(SCIPvarGetCol(subspacevars[j])), tableaurows[j], tableaurowinds[j], &ntableaurowinds[j]) );
269  	
270  	      if( ntableaurowinds[j] == -1 )
271  	      {
272  	         assert(sparse == 0 || sparse == -1);
273  	         sparse = 0;
274  	
275  	         for( i = nrows - 1; i >= 0; --i )
276  	            rownorm[i] += tableaurows[j][i] * tableaurows[j][i];
277  	      }
278  	      else
279  	      {
280  	         assert(sparse == 1 || sparse == -1);
281  	         sparse = 1;
282  	
283  	         /* allocate temporary memory */
284  	         if( usedrowids == NULL )
285  	         {
286  	            SCIP_CALL( SCIPallocBufferArray(scip, &rowids, nrows) );
287  	            SCIP_CALL( SCIPallocBufferArray(scip, &usedrowids, nrows) );
288  	            BMSclearMemoryArray(usedrowids, nrows);
289  	         }
290  	
291  	         for( i = ntableaurowinds[j] - 1; i >= 0; --i )
292  	         {
293  	            tableaurowind = tableaurowinds[j][i];
294  	            rownorm[tableaurowind] += tableaurows[j][tableaurowind] * tableaurows[j][tableaurowind];
295  	            assert(usedrowids != NULL);  /* for lint */
296  	            if( !usedrowids[tableaurowind] )
297  	            {
298  	               usedrowids[tableaurowind] = TRUE;
299  	               rowids[nrowids] = tableaurowind; /*lint !e644*/
300  	               ++nrowids;
301  	               assert(nrowids <= nrows);
302  	            }
303  	         }
304  	      }
305  	   }
306  	
307  	   /* compute ray direction (dense) */
308  	   if( sparse == 0 )
309  	   {
310  	      /* take average over all rows of the tableau */
311  	      for( i = nrows - 1; i >= 0; --i )
312  	      {
313  	         if( SCIPisFeasZero(scip, rownorm[i]) )
314  	            continue;
315  	         else
316  	            rownorm[i] = SQRT(rownorm[i]);
317  	
318  	         if( weighted )
319  	         {
320  	            rowweight = SCIProwGetDualsol(rows[i]);
321  	            if( SCIPisFeasZero(scip, rowweight) )
322  	               continue;
323  	         }
324  	         else
325  	            rowweight = 1.0;
326  	
327  	         for( j = nsubspacevars - 1; j >= 0; --j )
328  	         {
329  	            raydirection[j] += tableaurows[j][i] / (rownorm[i] * rowweight);
330  	            assert( ! SCIPisInfinity(scip, REALABS(raydirection[j])) );
331  	         }
332  	      }
333  	   }
334  	   /* compute ray direction (sparse) */
335  	   else
336  	   {
337  	      SCIP_Real* rowweights;
338  	      int r;
339  	      int k;
340  	
341  	      assert(usedrowids != NULL);
342  	      assert(rowids != NULL);
343  	      assert(sparse == 1);
344  	      assert(0 <= nrowids && nrowids <= nrows);
345  	
346  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowweights, nrows) );
347  	
348  	      /* compute norms of important rows and rowweights */
349  	      for( i = nrowids - 1; i >= 0; --i )
350  	      {
351  	         r = rowids[i];
352  	         assert(0 <= r && r < nrows);
353  	         assert(usedrowids[r]);
354  	
355  	         if( SCIPisFeasZero(scip, rownorm[r]) )
356  	         {
357  	            usedrowids[r] = FALSE;
358  	            --nrowids;
359  	            continue;
360  	         }
361  	         else
362  	            rownorm[r] = SQRT(rownorm[r]);
363  	
364  	         if( weighted )
365  	         {
366  	            rowweights[r] = SCIProwGetDualsol(rows[r]);
367  	            if( SCIPisFeasZero(scip, rowweights[r]) )
368  	            {
369  	               usedrowids[r] = FALSE;
370  	               --nrowids;
371  	               continue;
372  	            }
373  	         }
374  	         else
375  	            rowweights[r] = 1.0;
376  	      }
377  	
378  	      if( nrowids > 0 )
379  	      {
380  	         /* take average over all rows of the tableau */
381  	         for( j = nsubspacevars - 1; j >= 0; --j )
382  	         {
383  	            for( k = ntableaurowinds[j] - 1; k >= 0; --k )
384  	            {
385  	               tableaurowind = tableaurowinds[j][k];
386  	
387  	               if( usedrowids[tableaurowind] )
388  	               {
389  	                  raydirection[j] += tableaurows[j][tableaurowind] / (rownorm[tableaurowind] * rowweights[tableaurowind]);
390  	                  assert( ! SCIPisInfinity(scip, REALABS(raydirection[j])) );
391  	               }
392  	            }
393  	         }
394  	      }
395  	
396  	      SCIPfreeBufferArray(scip, &rowweights);
397  	      SCIPfreeBufferArray(scip, &usedrowids);
398  	      SCIPfreeBufferArray(scip, &rowids);
399  	   }
400  	   assert(usedrowids == NULL);
401  	
402  	   /* free memory */
403  	   SCIPfreeBufferArray(scip, &rownorm);
404  	   for( j = 0; j < nsubspacevars; ++j )
405  	   {
406  	      SCIPfreeBufferArray(scip, &tableaurowinds[j]);
407  	      SCIPfreeBufferArray(scip, &tableaurows[j]);
408  	   }
409  	   SCIPfreeBufferArray(scip, &ntableaurowinds);
410  	   SCIPfreeBufferArray(scip, &tableaurowinds);
411  	   SCIPfreeBufferArray(scip, &tableaurows);
412  	
413  	   return SCIP_OKAY;
414  	}
415  	
416  	
417  	/** generates the direction of the shooting ray as the average of the normalized non-basic vars and rows */
418  	static
419  	SCIP_RETCODE generateAverageNBRay(
420  	   SCIP*                 scip,               /**< SCIP data structure                   */
421  	   SCIP_Real*            raydirection,       /**< shooting ray                          */
422  	   int*                  fracspace,          /**< index set of fractional variables     */
423  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
424  	   int                   nsubspacevars       /**< dimension of fractional space         */
425  	   )
426  	{
427  	   SCIP_ROW** rows;
428  	   SCIP_COL** cols;
429  	   int nrows;
430  	   int ncols;
431  	   int i;
432  	
433  	   assert(scip != NULL);
434  	   assert(raydirection != NULL);
435  	   assert(fracspace != NULL);
436  	   assert(subspacevars != NULL);
437  	
438  	   SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
439  	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
440  	
441  	   /* add up non-basic variables */
442  	   for( i = nsubspacevars - 1; i >= 0; --i )
443  	   {
444  	      SCIP_Real solval;
445  	
446  	      solval = SCIPvarGetLPSol(subspacevars[i]);
447  	
448  	      if( SCIPisFeasEQ(scip, solval, SCIPvarGetLbLocal(subspacevars[i])) )
449  	         raydirection[i] = +1.0;
450  	      else if( SCIPisFeasEQ(scip, solval, SCIPvarGetUbLocal(subspacevars[i])) )
451  	         raydirection[i] = -1.0;
452  	      else
453  	         raydirection[i] = 0.0;
454  	   }
455  	
456  	   /* add up non-basic rows */
457  	   for( i = nrows - 1; i >= 0; --i )
458  	   {
459  	      SCIP_Real dualsol;
460  	      SCIP_Real factor;
461  	      SCIP_Real* coeffs;
462  	      SCIP_Real rownorm;
463  	      int j;
464  	      int nnonz;
465  	
466  	      dualsol = SCIProwGetDualsol(rows[i]);
467  	      if( SCIPisFeasPositive(scip, dualsol) )
468  	         factor = 1.0;
469  	      else if( SCIPisFeasNegative(scip, dualsol) )
470  	         factor = -1.0;
471  	      else
472  	         continue;
473  	
474  	      /* get the row's data */
475  	      coeffs = SCIProwGetVals(rows[i]);
476  	      cols = SCIProwGetCols(rows[i]);
477  	
478  	      nnonz = SCIProwGetNNonz(rows[i]);
479  	
480  	      rownorm = 0.0;
481  	      for( j = nnonz - 1; j >= 0; --j )
482  	      {
483  	         SCIP_VAR* var;
484  	         var = SCIPcolGetVar(cols[j]);
485  	         if( fracspace[SCIPvarGetProbindex(var)] >= 0 )
486  	            rownorm += coeffs[j] * coeffs[j];
487  	      }
488  	
489  	      if( SCIPisFeasZero(scip,rownorm) )
490  	         continue;
491  	      else
492  	      {
493  	         assert(rownorm > 0);
494  	         rownorm = SQRT(rownorm);
495  	      }
496  	
497  	      for( j = nnonz - 1; j >= 0; --j )
498  	      {
499  	         SCIP_VAR* var;
500  	         int f;
501  	
502  	         var = SCIPcolGetVar(cols[j]);
503  	         f = fracspace[SCIPvarGetProbindex(var)];
504  	
505  	         if( f >= 0 )
506  	         {
507  	            raydirection[f] += factor * coeffs[j] / rownorm;
508  	            assert( ! SCIPisInfinity(scip, REALABS(raydirection[f])) );
509  	         }
510  	      }
511  	   }
512  	   return SCIP_OKAY;
513  	}
514  	
515  	/** generates the starting point for the shooting ray in original coordinates */
516  	static
517  	void generateStartingPoint(
518  	   SCIP*                 scip,               /**< SCIP data structure                   */
519  	   SCIP_Real*            rayorigin,          /**< origin of the shooting ray            */
520  	   SCIP_VAR**            subspacevars,       /**< pointer to fractional space variables */
521  	   int                   nsubspacevars       /**< dimension of fractional space         */
522  	   )
523  	{
524  	   int v;
525  	
526  	   assert(scip != NULL);
527  	   assert(rayorigin != NULL);
528  	   assert(subspacevars != NULL);
529  	
530  	   for( v = nsubspacevars - 1; v >= 0; --v )
531  	      rayorigin[v] = SCIPvarGetLPSol(subspacevars[v]);
532  	}
533  	
534  	/** translates the inner point of the LP to an inner point rayorigin of the unit hyper octahedron and
535  	 *  transforms raydirection and rayorigin by reflections stored in sign
536  	 */
537  	static
538  	void flipCoords(
539  	   SCIP_Real*            rayorigin,          /**< origin of the shooting ray            */
540  	   SCIP_Real*            raydirection,       /**< direction of the shooting ray         */
541  	   SCIP_Bool*            sign,               /**< marker for flipped coordinates        */
542  	   int                   nsubspacevars       /**< dimension of fractional space         */
543  	   )
544  	{
545  	   int v;
546  	
547  	   assert(rayorigin != NULL);
548  	   assert(raydirection != NULL);
549  	   assert(sign != NULL);
550  	
551  	   for( v = nsubspacevars - 1; v >= 0; --v )
552  	   {
553  	      /* if raydirection[v] is negative, flip its sign */
554  	      if( raydirection[v] < 0 )
555  	      {
556  	         sign[v] = FALSE;
557  	         raydirection[v] *= -1.0;
558  	         rayorigin[v] *= -1.0; /* flip starting point in the same way like raydirection */
559  	      }
560  	      else
561  	         sign[v] = TRUE;
562  	   }
563  	}
564  	
565  	/** generates all facets, from which facet i could be obtained by a decreasing + to - flip
566  	 *  or a nonincreasing - to + flip and tests whether they are among the fmax nearest ones
567  	 */
568  	static
569  	void generateNeighborFacets(
570  	   SCIP*                 scip,               /**< SCIP data structure                   */
571  	   SCIP_Bool**           facets,             /**< facets got so far                     */
572  	   SCIP_Real*            lambda,             /**< distances of the facets               */
573  	   SCIP_Real*            rayorigin,          /**< origin of the shooting ray            */
574  	   SCIP_Real*            raydirection,       /**< direction of the shooting ray         */
575  	   SCIP_Real*            negquotient,        /**< array by which coordinates are sorted */
576  	   int                   nsubspacevars,      /**< dimension of fractional space         */
577  	   int                   f_max,              /**< maximal number of facets to create    */
578  	   int                   i,                  /**< current facet                         */
579  	   int*                  nfacets             /**< number of facets                      */
580  	   )
581  	{
582  	   SCIP_Real p;
583  	   SCIP_Real q;
584  	   SCIP_Real lam;
585  	   int minplus;
586  	   int j;
587  	
588  	   assert(scip != NULL);
589  	   assert(facets != NULL);
590  	   assert(facets[i] != NULL);
591  	   assert(lambda != NULL);
592  	   assert(rayorigin != NULL);
593  	   assert(raydirection != NULL);
594  	   assert(negquotient != NULL);
595  	   assert(nfacets != NULL);
596  	   assert(0 <= i && i < f_max);
597  	
598  	   /* determine the p and q values of the next facet to fix as a closest one */
599  	   p = 0.5 * nsubspacevars;
600  	   q = 0.0;
601  	   for( j = nsubspacevars - 1; j >= 0; --j )
602  	   {
603  	      if( facets[i][j] )
604  	      {
605  	         p -= rayorigin[j];
606  	         q += raydirection[j];
607  	      }
608  	      else
609  	      {
610  	         p += rayorigin[j];
611  	         q -= raydirection[j];
612  	      }
613  	   }
614  	
615  	   /* get the first + entry of the facet */
616  	   minplus = -1;
617  	   for( j = 0; j < nsubspacevars; ++j )
618  	   {
619  	      if( facets[i][j] )
620  	      {
621  	         minplus = j;
622  	         break;
623  	      }
624  	   }
625  	
626  	   /* facet (- - ... -) cannot be hit, because raydirection >= 0 */
627  	   assert(minplus >= 0);
628  	   assert(q != 0.0);
629  	   assert(SCIPisFeasEQ(scip, lambda[i], p/q));
630  	   assert(lambda[i] >= 0.0);
631  	
632  	   /* reverse search for facets from which the actual facet can be got by a single, decreasing + to - flip */
633  	   /* a facet will be inserted into the queue, iff it is one of the fmax closest ones already found */
634  	   for( j = 0; j < nsubspacevars && !facets[i][j] && SCIPisFeasGT(scip, negquotient[j], lambda[i]); ++j )
635  	   {
636  	      if( SCIPisFeasPositive(scip, q + 2*raydirection[j]) )
637  	      {
638  	         lam = (p - 2*rayorigin[j]) / (q + 2*raydirection[j]);
639  	         tryToInsert(scip, facets, lambda, i, j, f_max, nsubspacevars, lam, nfacets);
640  	      }
641  	   }
642  	
643  	   /* reverse search for facets from which the actual facet can be got by a single, nonincreasing - to + flip */
644  	   /* a facet will be inserted into the queue, iff it is one of the fmax closest ones already found */
645  	   for( j = nsubspacevars - 1; j >= 0 && facets[i][j] && SCIPisFeasLE(scip, negquotient[j], lambda[i]); --j )
646  	   {
647  	      if( SCIPisFeasPositive(scip, q - 2*raydirection[j]) )
648  	      {
649  	         lam = (p + 2*rayorigin[j]) / (q - 2*raydirection[j]);
650  	         if( negquotient[minplus] <= lam )
651  	            tryToInsert(scip, facets, lambda, i, j, f_max, nsubspacevars, lam, nfacets);
652  	      }
653  	   }
654  	#ifndef NDEBUG
655  	   for( j = 1; j < f_max; j++)
656  	      assert(SCIPisFeasGE(scip, lambda[j], lambda[j-1]));
657  	#endif
658  	}
659  	
660  	/** tests, whether an array is completely zero */
661  	static
662  	SCIP_Bool isZero(
663  	   SCIP*                 scip,               /**< SCIP data structure                   */
664  	   SCIP_Real*            raydirection,       /**< array to be checked                   */
665  	   int                   nsubspacevars       /**< size of array                         */
666  	   )
667  	{
668  	   int v;
669  	   SCIP_Bool iszero;
670  	
671  	   assert(scip != NULL);
672  	   assert(raydirection != NULL);
673  	   iszero = TRUE; 
674  	   for( v = nsubspacevars - 1; v >= 0; --v )
675  	   {
676  	      assert(!SCIPisInfinity(scip, raydirection[v]));
677  	
678  	      if( !SCIPisFeasZero(scip, raydirection[v]/100) )
679  	         iszero = FALSE;
680  	      else
681  	         raydirection[v] = 0.0;
682  	   }
683  	   return iszero;
684  	}
685  	
686  	
687  	/*
688  	 * Callback methods of primal heuristic
689  	 */
690  	
691  	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
692  	static
693  	SCIP_DECL_HEURCOPY(heurCopyOctane)
694  	{  /*lint --e{715}*/
695  	   assert(scip != NULL);
696  	   assert(heur != NULL);
697  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
698  	
699  	   /* call inclusion method of primal heuristic */
700  	   SCIP_CALL( SCIPincludeHeurOctane(scip) );
701  	
702  	   return SCIP_OKAY;
703  	}
704  	
705  	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
706  	static
707  	SCIP_DECL_HEURFREE(heurFreeOctane)
708  	{  /*lint --e{715}*/
709  	   SCIP_HEURDATA* heurdata;
710  	
711  	   assert(heur != NULL);
712  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
713  	   assert(scip != NULL);
714  	
715  	   /* free heuristic data */
716  	   heurdata = SCIPheurGetData(heur);
717  	   assert(heurdata != NULL);
718  	   SCIPfreeBlockMemory(scip, &heurdata);
719  	   SCIPheurSetData(heur, NULL);
720  	
721  	   return SCIP_OKAY;
722  	}
723  	
724  	
725  	/** initialization method of primal heuristic (called after problem was transformed) */
726  	static
727  	SCIP_DECL_HEURINIT(heurInitOctane)
728  	{  /*lint --e{715}*/
729  	   SCIP_HEURDATA* heurdata;
730  	
731  	   assert(heur != NULL);
732  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
733  	
734  	   /* get heuristic data */
735  	   heurdata = SCIPheurGetData(heur);
736  	   assert(heurdata != NULL);
737  	
738  	   /* create working solution */
739  	   SCIP_CALL( SCIPcreateSol(scip, &heurdata->sol, heur) );
740  	
741  	   /* initialize data */
742  	   heurdata->lastrule = 0;
743  	   heurdata->nsuccess = 0;
744  	
745  	   return SCIP_OKAY;
746  	}
747  	
748  	
749  	/** deinitialization method of primal heuristic (called before transformed problem is freed) */
750  	
751  	static
752  	SCIP_DECL_HEUREXIT(heurExitOctane)
753  	{  /*lint --e{715}*/
754  	   SCIP_HEURDATA* heurdata;
755  	
756  	   assert(heur != NULL);
757  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
758  	
759  	   /* get heuristic data */
760  	   heurdata = SCIPheurGetData(heur);
761  	   assert(heurdata != NULL);
762  	
763  	   /* free working solution */
764  	   SCIP_CALL( SCIPfreeSol(scip, &heurdata->sol) );
765  	
766  	   return SCIP_OKAY;
767  	}
768  	
769  	
770  	/** execution method of primal heuristic */
771  	static
772  	SCIP_DECL_HEUREXEC(heurExecOctane)
773  	{  /*lint --e{715}*/
774  	   SCIP_HEURDATA* heurdata;
775  	   SCIP_SOL* sol;
776  	   SCIP_SOL** first_sols;     /* stores the first ffirst sols in order to check for common violation of a row */
777  	
778  	   SCIP_VAR** vars;           /* the variables of the problem */
779  	   SCIP_VAR** fracvars;       /* variables, that are fractional in current LP solution */
780  	   SCIP_VAR** subspacevars;   /* the variables on which the search is performed. Either coinciding with vars or with the
781  	                               * space of all fractional variables of the current LP solution */
782  	
783  	   SCIP_Real p;               /* n/2 - <delta,x> ( for some facet delta ) */
784  	   SCIP_Real q;               /* <delta,a> */
785  	
786  	   SCIP_Real* rayorigin;      /* origin of the ray, vector x in paper */
787  	   SCIP_Real* raydirection;   /* direction of the ray, vector a in paper */
788  	   SCIP_Real* negquotient;    /* negated quotient of rayorigin and raydirection, vector v in paper */
789  	   SCIP_Real* lambda;         /* stores the distance of the facets (s.b.) to the origin of the ray */
790  	
791  	   SCIP_Bool usefracspace;    /* determines whether the search concentrates on fractional variables and fixes integer ones */
792  	   SCIP_Bool cons_viol;       /* used for checking whether a linear constraint is violated by one of the possible solutions */
793  	   SCIP_Bool success;
794  	   SCIP_Bool* sign;           /* signature of the direction of the ray */
795  	   SCIP_Bool** facets;        /* list of extended facets */
796  	
797  	   int nvars;            /* number of variables  */
798  	   int nbinvars;         /* number of 0-1-variables */
799  	   int nfracvars;        /* number of fractional variables in current LP solution */
800  	   int nsubspacevars;    /* dimension of the subspace on which the search is performed */
801  	   int nfacets;          /* number of facets hidden by the ray that where already found */
802  	   int i;                /* counter */
803  	   int j;                /* counter */
804  	   int f_max;            /* {0,1}-points to be checked */
805  	   int f_first;          /* {0,1}-points to be generated at first in order to check whether a restart is necessary */
806  	   int r;                /* counter */
807  	   int firstrule;
808  	
809  	   int* perm;            /* stores the way in which the coordinates were permuted */
810  	   int* fracspace;       /* maps the variables of the subspace to the original variables */
811  	
812  	   assert(heur != NULL);
813  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
814  	   assert(scip != NULL);
815  	   assert(result != NULL);
816  	   assert(SCIPhasCurrentNodeLP(scip));
817  	
818  	   *result = SCIP_DELAYED;
819  	
820  	   /* do not call heuristic of node was already detected to be infeasible */
821  	   if( nodeinfeasible )
822  	      return SCIP_OKAY;
823  	
824  	   /* only call heuristic, if an optimal LP solution is at hand */
825  	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
826  	      return SCIP_OKAY;
827  	
828  	   /* only call heuristic, if the LP objective value is smaller than the cutoff bound */
829  	   if( SCIPisGE(scip, SCIPgetLPObjval(scip), SCIPgetCutoffbound(scip)) )
830  	      return SCIP_OKAY;
831  	
832  	   *result = SCIP_DIDNOTRUN;
833  	
834  	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, &nbinvars, NULL, NULL, NULL) );
835  	
836  	   /* OCTANE is for use in 0-1 programs only */
837  	   if( nvars != nbinvars )
838  	      return SCIP_OKAY;
839  	
840  	   /* get heuristic's data */
841  	   heurdata = SCIPheurGetData(heur);
842  	   assert( heurdata != NULL );
843  	
844  	   /* don't call heuristic, if it was not successful enough in the past */
845  	   /*lint --e{647}*/
846  	   if( SCIPgetNNodes(scip) % (SCIPheurGetNCalls(heur) / (100 * SCIPheurGetNBestSolsFound(heur) + 10*heurdata->nsuccess + 1) + 1) != 0 )
847  	      return SCIP_OKAY;
848  	
849  	   SCIP_CALL( SCIPgetLPBranchCands(scip, &fracvars, NULL, NULL, &nfracvars, NULL, NULL) );
850  	
851  	   /* don't use integral starting points */
852  	   if( nfracvars == 0 )
853  	      return SCIP_OKAY;
854  	
855  	   /* get working pointers from heurdata */
856  	   sol = heurdata->sol;
857  	   assert( sol != NULL );
858  	   f_max = heurdata->f_max;
859  	   f_first = heurdata->f_first;
860  	   usefracspace = heurdata->usefracspace;
861  	
862  	   SCIP_CALL( SCIPallocBufferArray(scip, &fracspace, nvars) );
863  	
864  	   /* determine the space one which OCTANE should work either as the whole space or as the space of fractional variables */
865  	   if( usefracspace )
866  	   {
867  	      nsubspacevars = nfracvars;
868  	      SCIP_CALL( SCIPallocBufferArray(scip, &subspacevars, nsubspacevars) );
869  	      BMScopyMemoryArray(subspacevars, fracvars, nsubspacevars);
870  	      for( i = nvars - 1; i >= 0; --i )
871  	         fracspace[i] = -1;
872  	      for( i = nsubspacevars - 1; i >= 0; --i )
873  	         fracspace[SCIPvarGetProbindex(subspacevars[i])] = i;
874  	   }
875  	   else
876  	   {
877  	      int currentindex;
878  	
879  	      nsubspacevars = nvars;
880  	      SCIP_CALL( SCIPallocBufferArray(scip, &subspacevars, nsubspacevars) );
881  	
882  	      /* only copy the variables which are in the current LP */
883  	      currentindex = 0;
884  	      for( i = 0; i < nvars; ++i )
885  	      {
886  	         if( SCIPcolGetLPPos(SCIPvarGetCol(vars[i])) >= 0 )
887  	         {
888  	            subspacevars[currentindex] = vars[i];
889  	            fracspace[i] = currentindex;
890  	            ++currentindex;
891  	         }
892  	         else
893  	         {
894  	            fracspace[i] = -1;
895  	            --nsubspacevars;
896  	         }
897  	      }
898  	   }
899  	
900  	   /* nothing to do for empty search space */
901  	   if( nsubspacevars == 0 )
902  	   {
903  	      SCIPfreeBufferArray(scip, &subspacevars);
904  	      SCIPfreeBufferArray(scip, &fracspace);
905  	      return SCIP_OKAY;
906  	   }
907  	
908  	   assert(0 < nsubspacevars && nsubspacevars <= nvars);
909  	
910  	   for( i = 0; i < nsubspacevars; i++)
911  	      assert(fracspace[SCIPvarGetProbindex(subspacevars[i])] == i);
912  	
913  	   /* at most 2^(n-1) facets can be hit */
914  	   if( nsubspacevars < 30 )
915  	   {
916  	      /*lint --e{701}*/
917  	      assert(f_max > 0);
918  	      f_max = MIN(f_max, 1 << (nsubspacevars - 1) );
919  	   }
920  	
921  	   f_first = MIN(f_first, f_max);
922  	
923  	   /* memory allocation */
924  	   SCIP_CALL( SCIPallocBufferArray(scip, &rayorigin, nsubspacevars) );
925  	   SCIP_CALL( SCIPallocBufferArray(scip, &raydirection, nsubspacevars) );
926  	   SCIP_CALL( SCIPallocBufferArray(scip, &negquotient, nsubspacevars) );
927  	   SCIP_CALL( SCIPallocBufferArray(scip, &sign, nsubspacevars) );
928  	   SCIP_CALL( SCIPallocBufferArray(scip, &perm, nsubspacevars) );
929  	   SCIP_CALL( SCIPallocBufferArray(scip, &lambda, f_max + 1) );
930  	   SCIP_CALL( SCIPallocBufferArray(scip, &facets, f_max + 1) );
931  	
932  	   for( i = f_max; i >= 0; --i )
933  	   {
934  	      /*lint --e{866}*/
935  	      SCIP_CALL( SCIPallocBufferArray(scip, &facets[i], nsubspacevars) );
936  	   }
937  	   SCIP_CALL( SCIPallocBufferArray(scip, &first_sols, f_first) );
938  	
939  	   *result = SCIP_DIDNOTFIND;
940  	
941  	   /* starting OCTANE */
942  	   SCIPdebugMsg(scip, "run Octane heuristic on %s variables, which are %d vars, generate at most %d facets, using rule number %d\n",
943  	      usefracspace ? "fractional" : "all", nsubspacevars, f_max, (heurdata->lastrule+1)%5);
944  	
945  	   /* generate starting point in original coordinates */
946  	   generateStartingPoint(scip, rayorigin, subspacevars, nsubspacevars);
947  	   for( i = nsubspacevars - 1; i >= 0; --i )
948  	      rayorigin[i] -= 0.5;
949  	
950  	   firstrule = heurdata->lastrule;
951  	   ++firstrule;
952  	   for( r = firstrule; r <= firstrule + 5 && !SCIPisStopped(scip); r++ )
953  	   {
954  	      SCIP_ROW** rows;
955  	      int nrows;
956  	
957  	      /* generate shooting ray in original coordinates by certain rules */
958  	      switch(r % 5)
959  	      {
960  	      case 1:
961  	         if( !heurdata->useavgnbray )
962  	            continue;
963  	
964  	         SCIP_CALL( generateAverageNBRay(scip, raydirection, fracspace, subspacevars, nsubspacevars) );
965  	         break;
966  	      case 2:
967  	         if( !heurdata->useobjray )
968  	            continue;
969  	
970  	         SCIP_CALL( generateObjectiveRay(scip, raydirection, subspacevars, nsubspacevars) );
971  	         break;
972  	      case 3:
973  	         if( !heurdata->usediffray )
974  	            continue;
975  	
976  	         SCIP_CALL( generateDifferenceRay(scip, raydirection, subspacevars, nsubspacevars) );
977  	         break;
978  	      case 4:
979  	         if( !heurdata->useavgwgtray || !SCIPisLPSolBasic(scip) )
980  	            continue;
981  	
982  	         SCIP_CALL( generateAverageRay(scip, raydirection, subspacevars, nsubspacevars, TRUE) );
983  	         break;
984  	      case 0:
985  	         if( !heurdata->useavgray || !SCIPisLPSolBasic(scip) )
986  	            continue;
987  	
988  	         SCIP_CALL( generateAverageRay(scip, raydirection, subspacevars, nsubspacevars, FALSE) );
989  	         break;
990  	      default:
991  	         SCIPerrorMessage("invalid ray rule identifier\n");
992  	         SCIPABORT();
993  	      }
994  	
995  	      /* there must be a feasible direction for the shooting ray */
996  	      if( isZero(scip, raydirection, nsubspacevars) )
997  	         continue;
998  	
999  	      /* transform coordinates such that raydirection >= 0 */
1000 	      flipCoords(rayorigin, raydirection, sign, nsubspacevars);
1001 	
1002 	      for( i = f_max - 1; i >= 0; --i)
1003 	         lambda[i] = SCIPinfinity(scip);
1004 	
1005 	      /* calculate negquotient, initialize perm, facets[0], p, and q */
1006 	      p = 0.5 * nsubspacevars;
1007 	      q = 0.0;
1008 	      for( i = nsubspacevars - 1; i >= 0; --i )
1009 	      {
1010 	         /* calculate negquotient, the ratio of rayorigin and raydirection, paying special attention to the case raydirection[i] == 0 */
1011 	         if( SCIPisFeasZero(scip, raydirection[i]) )
1012 	         {
1013 	            if( rayorigin[i] < 0 )
1014 	               negquotient[i] = SCIPinfinity(scip);
1015 	            else
1016 	               negquotient[i] = -SCIPinfinity(scip);
1017 	         }
1018 	         else
1019 	            negquotient[i] = - (rayorigin[i] / raydirection[i]);
1020 	
1021 	         perm[i] = i;
1022 	
1023 	         /* initialization of facets[0] to the all-one facet with p and q its characteristic values */
1024 	         facets[0][i] = TRUE;
1025 	         p -= rayorigin[i];
1026 	         q += raydirection[i];
1027 	      }
1028 	
1029 	      assert(SCIPisPositive(scip, q));
1030 	
1031 	      /* resort the coordinates in nonincreasing order of negquotient */
1032 	      SCIPsortDownRealRealRealBoolPtr(negquotient, raydirection, rayorigin, sign, (void**) subspacevars, nsubspacevars);
1033 	
1034 	#ifndef NDEBUG
1035 	      for( i = 0; i < nsubspacevars; i++ )
1036 	      {
1037 	         assert( raydirection[i] >= 0 );
1038 	         assert(!SCIPisInfinity(scip, REALABS(raydirection[i])));
1039 	      }
1040 	      for( i = 1; i < nsubspacevars; i++ )
1041 	         assert( negquotient[i - 1] >= negquotient[i] );
1042 	#endif
1043 	      /* finished initialization */
1044 	
1045 	      /* find the first facet of the octahedron hit by a ray shot from rayorigin into direction raydirection */
1046 	      for( i = 0; i < nsubspacevars && negquotient[i] * q > p; ++i )
1047 	      {
1048 	         facets[0][i] = FALSE;
1049 	         p += 2 * rayorigin[i];
1050 	         q -= 2 * raydirection[i];
1051 	         assert(SCIPisPositive(scip, p));
1052 	         assert(SCIPisPositive(scip, q));
1053 	      }
1054 	
1055 	      /* avoid dividing by values close to 0.0 */
1056 	      if( !SCIPisFeasPositive(scip, q) )
1057 	         continue;
1058 	
1059 	      /* assert necessary for flexelint */
1060 	      assert(q != 0.0);
1061 	      lambda[0] = p / q;
1062 	
1063 	      nfacets = 1;
1064 	
1065 	      /* find the first facets hit by the ray */
1066 	      for( i = 0; i < nfacets && i < f_first; ++i)
1067 	         generateNeighborFacets(scip, facets, lambda, rayorigin, raydirection, negquotient, nsubspacevars, f_max, i, &nfacets);
1068 	
1069 	      /* construct the first ffirst possible solutions */
1070 	      for( i = 0; i < nfacets && i < f_first; ++i )
1071 	      {
1072 	         SCIP_CALL( SCIPcreateSol(scip, &first_sols[i], heur) );
1073 	         SCIP_CALL( getSolFromFacet(scip, facets[i], first_sols[i], sign, subspacevars, nsubspacevars) );
1074 	         assert( first_sols[i] != NULL );
1075 	      }
1076 	
1077 	      /* try, whether there is a row violated by all of the first ffirst solutions */
1078 	      cons_viol = FALSE;
1079 	      SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
1080 	      for( i = nrows - 1; i >= 0; --i )
1081 	      {
1082 	         if( !SCIProwIsLocal(rows[i]) )
1083 	         {
1084 	            SCIP_COL** cols;
1085 	            SCIP_Real constant;
1086 	            SCIP_Real lhs;
1087 	            SCIP_Real rhs;
1088 	            SCIP_Real rowval;
1089 	            SCIP_Real* coeffs;
1090 	            int nnonzerovars;
1091 	            int k;
1092 	
1093 	            /* get the row's data */
1094 	            constant = SCIProwGetConstant(rows[i]);
1095 	            lhs = SCIProwGetLhs(rows[i]);
1096 	            rhs = SCIProwGetRhs(rows[i]);
1097 	            coeffs = SCIProwGetVals(rows[i]);
1098 	            nnonzerovars = SCIProwGetNNonz(rows[i]);
1099 	            cols = SCIProwGetCols(rows[i]);
1100 	            rowval = constant;
1101 	
1102 	            for( j = nnonzerovars - 1; j >= 0; --j )
1103 	               rowval += coeffs[j] * SCIPgetSolVal(scip, first_sols[0], SCIPcolGetVar(cols[j]));
1104 	
1105 	            /* if the row's lhs is violated by the first sol, test, whether it is violated by the next ones, too */
1106 	            if( lhs > rowval )
1107 	            {
1108 	               cons_viol = TRUE;
1109 	               for( k = MIN(f_first, nfacets) - 1; k > 0; --k )
1110 	               {
1111 	                  rowval = constant;
1112 	                  for( j = nnonzerovars - 1; j >= 0; --j )
1113 	                     rowval += coeffs[j] * SCIPgetSolVal(scip, first_sols[k], SCIPcolGetVar(cols[j]));
1114 	                  if( lhs <= rowval )
1115 	                  {
1116 	                     cons_viol = FALSE;
1117 	                     break;
1118 	                  }
1119 	               }
1120 	            }
1121 	            /* dito for the right hand side */
1122 	            else if( rhs < rowval )
1123 	            {
1124 	               cons_viol = TRUE;
1125 	               for( k = MIN(f_first, nfacets) - 1; k > 0; --k )
1126 	               {
1127 	                  rowval = constant;
1128 	                  for( j = nnonzerovars - 1; j >= 0; --j )
1129 	                     rowval += coeffs[j] * SCIPgetSolVal(scip, first_sols[k], SCIPcolGetVar(cols[j]));
1130 	                  if( rhs >= rowval )
1131 	                  {
1132 	                     cons_viol = FALSE;
1133 	                     break;
1134 	                  }
1135 	               }
1136 	            }
1137 	            /* break as soon as one row is violated by all of the ffirst solutions */
1138 	            if( cons_viol )
1139 	               break;
1140 	         }
1141 	      }
1142 	
1143 	      if( !cons_viol )
1144 	      {
1145 	         /* if there was no row violated by all solutions, try whether one or more of them are feasible */
1146 	         for( i = MIN(f_first, nfacets) - 1; i >= 0; --i )
1147 	         {
1148 	            assert(first_sols[i] != NULL);
1149 	            SCIP_CALL( SCIPtrySol(scip, first_sols[i], FALSE, FALSE, TRUE, FALSE, TRUE, &success) );
1150 	            if( success )
1151 	               *result = SCIP_FOUNDSOL;
1152 	         }
1153 	         /* search for further facets and construct and try solutions out of facets fixed as closest ones */
1154 	         for( i = f_first; i < f_max; ++i)
1155 	         {
1156 	            if( i >= nfacets )
1157 	               break;
1158 	            generateNeighborFacets(scip, facets, lambda, rayorigin, raydirection, negquotient, nsubspacevars, f_max, i, &nfacets);
1159 	            SCIP_CALL( getSolFromFacet(scip, facets[i], sol, sign, subspacevars, nsubspacevars) );
1160 	            SCIP_CALL( SCIPtrySol(scip, sol, FALSE, FALSE, TRUE, FALSE, TRUE, &success) );
1161 	            if( success )
1162 	               *result = SCIP_FOUNDSOL;
1163 	         }
1164 	      }
1165 	
1166 	      /* finished OCTANE */
1167 	      for( i = MIN(f_first, nfacets) - 1; i >= 0; --i )
1168 	      {
1169 	         SCIP_CALL( SCIPfreeSol(scip, &first_sols[i]) );
1170 	      }
1171 	   }
1172 	   heurdata->lastrule = r;
1173 	
1174 	   if( *result == SCIP_FOUNDSOL )
1175 	      ++(heurdata->nsuccess);
1176 	
1177 	   /* free temporary memory */
1178 	   SCIPfreeBufferArray(scip, &first_sols);
1179 	   for( i = 0; i <= f_max; ++i )
1180 	      SCIPfreeBufferArray(scip, &facets[i]);
1181 	   SCIPfreeBufferArray(scip, &facets);
1182 	   SCIPfreeBufferArray(scip, &lambda);
1183 	   SCIPfreeBufferArray(scip, &perm);
1184 	   SCIPfreeBufferArray(scip, &sign);
1185 	   SCIPfreeBufferArray(scip, &negquotient);
1186 	   SCIPfreeBufferArray(scip, &raydirection);
1187 	   SCIPfreeBufferArray(scip, &rayorigin);
1188 	   SCIPfreeBufferArray(scip, &subspacevars);
1189 	   SCIPfreeBufferArray(scip, &fracspace);
1190 	
1191 	   return SCIP_OKAY;
1192 	}
1193 	
1194 	
1195 	/*
1196 	 * primal heuristic specific interface methods
1197 	 */
1198 	
1199 	/** creates the octane primal heuristic and includes it in SCIP */
1200 	SCIP_RETCODE SCIPincludeHeurOctane(
1201 	   SCIP*                 scip                /**< SCIP data structure */
1202 	   )
1203 	{
1204 	   SCIP_HEURDATA* heurdata;
1205 	   SCIP_HEUR* heur;
1206 	
1207 	   /* create Octane primal heuristic data */
1208 	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
1209 	
1210 	   /* include primal heuristic */
1211 	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
1212 	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
1213 	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecOctane, heurdata) );
1214 	
1215 	   assert(heur != NULL);
1216 	
1217 	   /* set non-NULL pointers to callback methods */
1218 	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyOctane) );
1219 	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeOctane) );
1220 	   SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitOctane) );
1221 	   SCIP_CALL( SCIPsetHeurExit(scip, heur, heurExitOctane) );
1222 	
1223 	   /* add octane primal heuristic parameters */
1224 	   SCIP_CALL( SCIPaddIntParam(scip,
1225 	         "heuristics/octane/fmax",
1226 	         "number of 0-1-points to be tested as possible solutions by OCTANE",
1227 	         &heurdata->f_max, TRUE, DEFAULT_FMAX, 1, INT_MAX, NULL, NULL) );
1228 	
1229 	   SCIP_CALL( SCIPaddIntParam(scip,
1230 	         "heuristics/octane/ffirst",
1231 	         "number of 0-1-points to be tested at first whether they violate a common row",
1232 	         &heurdata->f_first, TRUE, DEFAULT_FFIRST, 1, INT_MAX, NULL, NULL) );
1233 	
1234 	   SCIP_CALL( SCIPaddBoolParam(scip,
1235 	         "heuristics/octane/usefracspace",
1236 	         "execute OCTANE only in the space of fractional variables (TRUE) or in the full space?",
1237 	         &heurdata->usefracspace, TRUE, DEFAULT_USEFRACSPACE, NULL, NULL) );
1238 	
1239 	   SCIP_CALL( SCIPaddBoolParam(scip,
1240 	         "heuristics/octane/useobjray",
1241 	         "should the inner normal of the objective be used as one ray direction?",
1242 	         &heurdata->useobjray, TRUE, TRUE, NULL, NULL) );
1243 	
1244 	   SCIP_CALL( SCIPaddBoolParam(scip,
1245 	         "heuristics/octane/useavgray",
1246 	         "should the average of the basic cone be used as one ray direction?",
1247 	         &heurdata->useavgray, TRUE, TRUE, NULL, NULL) );
1248 	
1249 	   SCIP_CALL( SCIPaddBoolParam(scip,
1250 	         "heuristics/octane/usediffray",
1251 	         "should the difference between the root solution and the current LP solution be used as one ray direction?",
1252 	         &heurdata->usediffray, TRUE, FALSE, NULL, NULL) );
1253 	
1254 	   SCIP_CALL( SCIPaddBoolParam(scip,
1255 	         "heuristics/octane/useavgwgtray",
1256 	         "should the weighted average of the basic cone be used as one ray direction?",
1257 	         &heurdata->useavgwgtray, TRUE, TRUE, NULL, NULL) );
1258 	
1259 	   SCIP_CALL( SCIPaddBoolParam(scip,
1260 	         "heuristics/octane/useavgnbray",
1261 	         "should the weighted average of the nonbasic cone be used as one ray direction?",
1262 	         &heurdata->useavgnbray, TRUE, TRUE, NULL, NULL) );
1263 	
1264 	   return SCIP_OKAY;
1265 	}
1266