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_ofins.c
26   	 * @ingroup DEFPLUGINS_HEUR
27   	 * @brief  OFINS - Objective Function Induced Neighborhood Search - a primal heuristic for reoptimization
28   	 * @author Jakob Witzig
29   	 */
30   	
31   	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
32   	
33   	#include "blockmemshell/memory.h"
34   	#include "scip/heuristics.h"
35   	#include "scip/heur_ofins.h"
36   	#include "scip/pub_event.h"
37   	#include "scip/pub_heur.h"
38   	#include "scip/pub_message.h"
39   	#include "scip/pub_misc.h"
40   	#include "scip/pub_sol.h"
41   	#include "scip/pub_var.h"
42   	#include "scip/scip_branch.h"
43   	#include "scip/scip_copy.h"
44   	#include "scip/scip_event.h"
45   	#include "scip/scip_general.h"
46   	#include "scip/scip_heur.h"
47   	#include "scip/scip_mem.h"
48   	#include "scip/scip_message.h"
49   	#include "scip/scip_nodesel.h"
50   	#include "scip/scip_numerics.h"
51   	#include "scip/scip_param.h"
52   	#include "scip/scip_prob.h"
53   	#include "scip/scip_sol.h"
54   	#include "scip/scip_solve.h"
55   	#include "scip/scip_solvingstats.h"
56   	#include "scip/scip_timing.h"
57   	#include <string.h>
58   	
59   	#define HEUR_NAME             "ofins"
60   	#define HEUR_DESC             "primal heuristic for reoptimization, objective function induced neighborhood search"
61   	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_LNS
62   	#define HEUR_PRIORITY         60000
63   	#define HEUR_FREQ             0
64   	#define HEUR_FREQOFS          0
65   	#define HEUR_MAXDEPTH         0
66   	#define HEUR_TIMING           SCIP_HEURTIMING_BEFORENODE
67   	#define HEUR_USESSUBSCIP      TRUE  /**< does the heuristic use a secondary SCIP instance? */
68   	
69   	/* default values for OFINS-specific plugins */
70   	#define DEFAULT_MAXNODES      5000LL    /**< maximum number of nodes to regard in the subproblem */
71   	#define DEFAULT_MAXCHGRATE    0.50      /**< maximum percentage of changed objective coefficients */
72   	#define DEFAULT_COPYCUTS      TRUE      /**< if DEFAULT_USELPROWS is FALSE, then should all active cuts from the cutpool
73   	                                         *   of the original scip be copied to constraints of the subscip */
74   	#define DEFAULT_MAXCHANGE     0.04      /**< maximal rate of change per coefficient to get fixed */
75   	#define DEFAULT_MINIMPROVE    0.01      /**< factor by which OFINS should at least improve the incumbent */
76   	#define DEFAULT_ADDALLSOLS   FALSE      /**< should all subproblem solutions be added to the original SCIP? */
77   	#define DEFAULT_MINNODES      50LL      /**< minimum number of nodes to regard in the subproblem */
78   	#define DEFAULT_NODESOFS      500LL     /**< number of nodes added to the contingent of the total nodes */
79   	#define DEFAULT_NODESQUOT     0.1       /**< subproblem nodes in relation to nodes of the original problem */
80   	#define DEFAULT_LPLIMFAC      2.0       /**< factor by which the limit on the number of LP depends on the node limit */
81   	
82   	/* event handler properties */
83   	#define EVENTHDLR_NAME         "Ofins"
84   	#define EVENTHDLR_DESC         "LP event handler for " HEUR_NAME " heuristic"
85   	
86   	
87   	/** primal heuristic data */
88   	struct SCIP_HeurData
89   	{
90   	   SCIP_Real             maxchangerate;      /**< maximal rate of changed coefficients in the objective function */
91   	   SCIP_Longint          maxnodes;           /**< maximum number of nodes to regard in the subproblem */
92   	   SCIP_Bool             copycuts;           /**< should all active cuts from cutpool be copied to constraints in subproblem? */
93   	   SCIP_Bool             addallsols;         /**< should all subproblem solutions be added to the original SCIP? */
94   	   SCIP_Longint          minnodes;           /**< minimum number of nodes to regard in the subproblem */
95   	   SCIP_Longint          nodesofs;           /**< number of nodes added to the contingent of the total nodes */
96   	   SCIP_Real             maxchange;          /**< maximal rate of change per coefficient to get fixed */
97   	   SCIP_Real             minimprove;         /**< factor by which OFINS should at least improve the incumbent */
98   	   SCIP_Real             nodesquot;          /**< subproblem nodes in relation to nodes of the original problem */
99   	   SCIP_Real             nodelimit;          /**< the nodelimit employed in the current sub-SCIP, for the event handler*/
100  	   SCIP_Real             lplimfac;           /**< factor by which the limit on the number of LP depends on the node limit */
101  	};
102  	
103  	/* ---------------- Callback methods of event handler ---------------- */
104  	
105  	/* exec the event handler
106  	 *
107  	 * we interrupt the solution process
108  	 */
109  	static
110  	SCIP_DECL_EVENTEXEC(eventExecOfins)
111  	{
112  	   SCIP_HEURDATA* heurdata;
113  	
114  	   assert(eventhdlr != NULL);
115  	   assert(eventdata != NULL);
116  	   assert(strcmp(SCIPeventhdlrGetName(eventhdlr), EVENTHDLR_NAME) == 0);
117  	   assert(event != NULL);
118  	   assert(SCIPeventGetType(event) & SCIP_EVENTTYPE_LPSOLVED);
119  	
120  	   heurdata = (SCIP_HEURDATA*)eventdata;
121  	   assert(heurdata != NULL);
122  	
123  	   /* interrupt solution process of sub-SCIP */
124  	   if( SCIPgetNLPs(scip) > heurdata->lplimfac * heurdata->nodelimit )
125  	   {
126  	      SCIPdebugMsg(scip, "interrupt after %" SCIP_LONGINT_FORMAT " LPs\n",SCIPgetNLPs(scip));
127  	      SCIP_CALL( SCIPinterruptSolve(scip) );
128  	   }
129  	
130  	   return SCIP_OKAY;
131  	}
132  	
133  	/* setup and solve the sub-SCIP */
134  	static
135  	SCIP_RETCODE setupAndSolve(
136  	   SCIP*                 scip,               /**< original SCIP data structure */
137  	   SCIP*                 subscip,            /**< sub-SCIP data structure */
138  	   SCIP_HEUR*            heur,               /**< heuristic data structure */
139  	   SCIP_HEURDATA*        heurdata,           /**< euristic's private data structure */
140  	   SCIP_RESULT*          result,             /**< result data structure */
141  	   SCIP_Longint          nstallnodes,        /**< number of stalling nodes for the subproblem */
142  	   SCIP_Bool*            chgcoeffs           /**< array of changed coefficients */
143  	
144  	   )
145  	{
146  	   SCIP_HASHMAP* varmapfw;
147  	   SCIP_VAR** vars;
148  	   SCIP_VAR** subvars;
149  	   SCIP_EVENTHDLR* eventhdlr;
150  	
151  	   SCIP_SOL* sol;
152  	   SCIP_VAR** fixedvars;
153  	   SCIP_Real* fixedvals;
154  	   int nfixedvars;
155  	
156  	   int nvars;
157  	   int nintvars;
158  	   int i;
159  	
160  	   SCIP_SOL** subsols;
161  	   int nsubsols = 0;
162  	
163  	   SCIP_Bool success;
164  	   SCIP_RETCODE retcode;
165  	   SCIP_STATUS status;
166  	
167  	   assert(scip != NULL);
168  	   assert(subscip != NULL);
169  	   assert(heur != NULL);
170  	   assert(heurdata != NULL);
171  	   assert(result != NULL);
172  	   assert(chgcoeffs != NULL);
173  	
174  	   SCIPdebugMsg(scip, "+---+ Start OFINS heuristic +---+\n");
175  	
176  	   /* get variable data */
177  	   vars = SCIPgetVars(scip);
178  	   nvars = SCIPgetNVars(scip);
179  	
180  	   /* create the variable mapping hash map */
181  	   SCIP_CALL( SCIPhashmapCreate(&varmapfw, SCIPblkmem(subscip), nvars) );
182  	
183  	   /* get optimal solution of the last iteration */
184  	   sol = SCIPgetReoptLastOptSol(scip);
185  	
186  	   /* if the solution is NULL the last problem was infeasible */
187  	   if( sol == NULL )
188  	      return SCIP_OKAY;
189  	
190  	   nintvars = SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) + SCIPgetNImplVars(scip);
191  	   SCIP_CALL( SCIPallocBufferArray(scip, &fixedvars, nvars) );
192  	   SCIP_CALL( SCIPallocBufferArray(scip, &fixedvals, nvars) );
193  	
194  	   /* determine variables to fix in the sub-SCIP */
195  	   nfixedvars = 0;
196  	   for( i = 0; i < nintvars; i++ )
197  	   {
198  	      if( !chgcoeffs[i] )
199  	      {
200  	         fixedvars[nfixedvars] = vars[i];
201  	         fixedvals[nfixedvars] = SCIPgetSolVal(scip, sol, vars[i]);
202  	         ++nfixedvars;
203  	      }
204  	   }
205  	
206  	   /* create a problem copy as sub SCIP */
207  	   SCIP_CALL( SCIPcopyLargeNeighborhoodSearch(scip, subscip, varmapfw, "ofins", fixedvars, fixedvals, nfixedvars, FALSE,
208  	         FALSE, &success, NULL) );
209  	   assert(success);
210  	
211  	   SCIPfreeBufferArrayNull(scip, &fixedvals);
212  	   SCIPfreeBufferArrayNull(scip, &fixedvars);
213  	
214  	   /* create event handler for LP events */
215  	   eventhdlr = NULL;
216  	   SCIP_CALL( SCIPincludeEventhdlrBasic(subscip, &eventhdlr, EVENTHDLR_NAME, EVENTHDLR_DESC, eventExecOfins, NULL) );
217  	   if( eventhdlr == NULL )
218  	   {
219  	      SCIPerrorMessage("event handler for " HEUR_NAME " heuristic not found.\n");
220  	      return SCIP_PLUGINNOTFOUND;
221  	   }
222  	
223  	   SCIP_CALL( SCIPallocBufferArray(scip, &subvars, nvars) );
224  	   for( i = 0; i < nvars; i++ )
225  	     subvars[i] = (SCIP_VAR*) SCIPhashmapGetImage(varmapfw, vars[i]);
226  	
227  	   /* free hash map */
228  	   SCIPhashmapFree(&varmapfw);
229  	
230  	   /* set an objective limit */
231  	   SCIPdebugMsg(scip, "set objective limit of %g to sub-SCIP\n", SCIPgetUpperbound(scip));
232  	   SCIP_CALL( SCIPsetObjlimit(subscip, SCIPgetUpperbound(scip)) );
233  	
234  	   SCIPdebugMsg(scip, "OFINS subproblem: %d vars, %d cons\n", SCIPgetNVars(subscip), SCIPgetNConss(subscip));
235  	
236  	   /* do not abort subproblem on CTRL-C */
237  	   SCIP_CALL( SCIPsetBoolParam(subscip, "misc/catchctrlc", FALSE) );
238  	
239  	#ifdef SCIP_DEBUG
240  	   /* for debugging, enable full output */
241  	   SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 5) );
242  	   SCIP_CALL( SCIPsetIntParam(subscip, "display/freq", 100000000) );
243  	#else
244  	   /* disable statistic timing inside sub SCIP and output to console */
245  	   SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
246  	   SCIP_CALL( SCIPsetBoolParam(subscip, "timing/statistictiming", FALSE) );
247  	#endif
248  	
249  	   /* set limits for the subproblem */
250  	   SCIP_CALL( SCIPcopyLimits(scip, subscip) );
251  	   heurdata->nodelimit = heurdata->maxnodes;
252  	   SCIP_CALL( SCIPsetLongintParam(subscip, "limits/stallnodes", nstallnodes) );
253  	   SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", heurdata->maxnodes) );
254  	
255  	   /* forbid recursive call of heuristics and separators solving sub-SCIPs */
256  	   SCIP_CALL( SCIPsetSubscipsOff(subscip, TRUE) );
257  	
258  	   /* disable cutting plane separation */
259  	   SCIP_CALL( SCIPsetSeparating(subscip, SCIP_PARAMSETTING_OFF, TRUE) );
260  	
261  	   /* disable expensive presolving */
262  	   SCIP_CALL( SCIPsetPresolving(subscip, SCIP_PARAMSETTING_FAST, TRUE) );
263  	
264  	   /* use best estimate node selection */
265  	   if( SCIPfindNodesel(subscip, "estimate") != NULL && !SCIPisParamFixed(subscip, "nodeselection/estimate/stdpriority") )
266  	   {
267  	      SCIP_CALL( SCIPsetIntParam(subscip, "nodeselection/estimate/stdpriority", INT_MAX/4) );
268  	   }
269  	
270  	   /* use inference branching */
271  	   if( SCIPfindBranchrule(subscip, "inference") != NULL && !SCIPisParamFixed(subscip, "branching/inference/priority") )
272  	   {
273  	      SCIP_CALL( SCIPsetIntParam(subscip, "branching/inference/priority", INT_MAX/4) );
274  	   }
275  	
276  	   /* disable conflict analysis */
277  	   if( !SCIPisParamFixed(subscip, "conflict/enable") )
278  	   {
279  	      SCIP_CALL( SCIPsetBoolParam(subscip, "conflict/enable", FALSE) );
280  	   }
281  	
282  	   /* speed up sub-SCIP by not checking dual LP feasibility */
283  	   SCIP_CALL( SCIPsetBoolParam(subscip, "lp/checkdualfeas", FALSE) );
284  	
285  	   /* presolve the subproblem */
286  	   retcode = SCIPpresolve(subscip);
287  	
288  	   /* errors in solving the subproblem should not kill the overall solving process;
289  	    * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop.
290  	    */
291  	   if( retcode != SCIP_OKAY )
292  	   {
293  	      SCIPwarningMessage(scip, "Error while presolving subproblem in %s heuristic; sub-SCIP terminated with code <%d>\n", HEUR_NAME, retcode);
294  	
295  	      SCIPABORT(); /*lint --e{527}*/
296  	
297  	      /* free */
298  	      SCIPfreeBufferArray(scip, &subvars);
299  	      return SCIP_OKAY;
300  	   }
301  	
302  	   SCIPdebugMsg(scip, "%s presolved subproblem: %d vars, %d cons\n", HEUR_NAME, SCIPgetNVars(subscip), SCIPgetNConss(subscip));
303  	
304  	   assert(eventhdlr != NULL);
305  	
306  	   SCIP_CALL( SCIPtransformProb(subscip) );
307  	   SCIP_CALL( SCIPcatchEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, NULL) );
308  	
309  	   /* solve the subproblem */
310  	   SCIPdebugMsg(scip, "solving subproblem: nstallnodes=%" SCIP_LONGINT_FORMAT ", maxnodes=%" SCIP_LONGINT_FORMAT "\n", nstallnodes, heurdata->maxnodes);
311  	
312  	   /* errors in solving the subproblem should not kill the overall solving process;
313  	    * hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop.
314  	    */
315  	   SCIP_CALL_ABORT( SCIPsolve(subscip) );
316  	
317  	   SCIP_CALL( SCIPdropEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, -1) );
318  	
319  	   /* print solving statistics of subproblem if we are in SCIP's debug mode */
320  	   SCIPdebug( SCIP_CALL( SCIPprintStatistics(subscip, NULL) ) );
321  	
322  	   status = SCIPgetStatus(subscip);
323  	
324  	   switch (status) {
325  	   case SCIP_STATUS_INFEASIBLE:
326  	      break;
327  	   case SCIP_STATUS_INFORUNBD:
328  	   case SCIP_STATUS_UNBOUNDED:
329  	   {
330  	      /* transfer the primal ray from the sub-SCIP to the main SCIP */
331  	      if( SCIPhasPrimalRay(subscip) )
332  	      {
333  	         SCIP_SOL* primalray;
334  	
335  	         SCIP_CALL( SCIPcreateSol(scip, &primalray, heur) );
336  	
337  	         /* transform the ray into the space of the source scip */
338  	         for( i = 0; i < nvars; i++ )
339  	         {
340  	            SCIP_CALL( SCIPsetSolVal(scip, primalray, vars[i],
341  	                  subvars[i] != NULL ? SCIPgetPrimalRayVal(subscip, subvars[i]) : 0.0) );
342  	         }
343  	
344  	         SCIPdebug( SCIP_CALL( SCIPprintRay(scip, primalray, 0, FALSE) ); );
345  	
346  	         /* update the primal ray of the source scip */
347  	         SCIP_CALL( SCIPupdatePrimalRay(scip, primalray) );
348  	         SCIP_CALL( SCIPfreeSol(scip, &primalray) );
349  	
350  	         *result = SCIP_UNBOUNDED;
351  	      }
352  	
353  	      break;
354  	   }
355  	   default:
356  	      /* check, whether a solution was found;
357  	       * due to numerics, it might happen that not all solutions are feasible -> try all solutions until one was accepted
358  	       */
359  	      nsubsols = SCIPgetNSols(subscip);
360  	      subsols = SCIPgetSols(subscip);
361  	      success = FALSE;
362  	      for( i = 0; i < nsubsols && (!success || heurdata->addallsols); i++ )
363  	      {
364  	         SCIP_SOL* newsol;
365  	
366  	         SCIP_CALL( SCIPtranslateSubSol(scip, subscip, subsols[i], heur, subvars, &newsol) );
367  	
368  	         /* try to add new solution to scip and free it immediately */
369  	         SCIP_CALL( SCIPtrySolFree(scip, &newsol, FALSE, FALSE, TRUE, TRUE, TRUE, &success) );
370  	
371  	         if( success )
372  	            *result = SCIP_FOUNDSOL;
373  	      }
374  	      break;
375  	   } /*lint !e788*/
376  	
377  	   SCIPstatisticPrintf("%s statistic: fixed %6.3f integer variables, needed %6.1f seconds, %" SCIP_LONGINT_FORMAT " nodes, solution %10.4f found at node %" SCIP_LONGINT_FORMAT "\n",
378  	      HEUR_NAME, 0.0, SCIPgetSolvingTime(subscip), SCIPgetNNodes(subscip), success ? SCIPgetPrimalbound(scip) : SCIPinfinity(scip),
379  	      nsubsols > 0 ? SCIPsolGetNodenum(SCIPgetBestSol(subscip)) : -1 );
380  	
381  	   /* free subproblem */
382  	   SCIPfreeBufferArray(scip, &subvars);
383  	
384  	   return SCIP_OKAY;
385  	}
386  	
387  	/** main procedure of the OFINS heuristic, creates and solves a sub-SCIP */
388  	static
389  	SCIP_RETCODE applyOfins(
390  	   SCIP*                 scip,               /**< original SCIP data structure */
391  	   SCIP_HEUR*            heur,               /**< heuristic data structure */
392  	   SCIP_HEURDATA*        heurdata,           /**< euristic's private data structure */
393  	   SCIP_RESULT*          result,             /**< result data structure */
394  	   SCIP_Longint          nstallnodes,        /**< number of stalling nodes for the subproblem */
395  	   SCIP_Bool*            chgcoeffs           /**< array of changed coefficients */
396  	   )
397  	{
398  	   SCIP* subscip;
399  	   SCIP_RETCODE retcode;
400  	   SCIP_Bool success;
401  	
402  	   assert(scip != NULL);
403  	   assert(heur != NULL);
404  	   assert(heurdata != NULL);
405  	   assert(result != NULL);
406  	   assert(chgcoeffs != NULL);
407  	
408  	   *result = SCIP_DIDNOTRUN;
409  	
410  	   /* check whether there is enough time and memory left */
411  	   SCIP_CALL( SCIPcheckCopyLimits(scip, &success) );
412  	
413  	   if( !success )
414  	      return SCIP_OKAY;
415  	
416  	   *result = SCIP_DIDNOTFIND;
417  	
418  	   /* do not run, if no solution was found */
419  	   if ( SCIPgetReoptLastOptSol(scip) == NULL )
420  	      return SCIP_OKAY;
421  	
422  	   /* initialize the subproblem */
423  	   SCIP_CALL( SCIPcreate(&subscip) );
424  	
425  	   retcode = setupAndSolve(scip, subscip, heur, heurdata, result, nstallnodes, chgcoeffs);
426  	
427  	   SCIP_CALL( SCIPfree(&subscip) );
428  	
429  	   SCIP_CALL( retcode );
430  	
431  	   return SCIP_OKAY;
432  	}
433  	
434  	
435  	
436  	
437  	/*
438  	 * Callback methods of primal heuristic
439  	 */
440  	
441  	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
442  	static
443  	SCIP_DECL_HEURCOPY(heurCopyOfins)
444  	{  /*lint --e{715}*/
445  	   assert(scip != NULL);
446  	   assert(heur != NULL);
447  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
448  	
449  	   /* call inclusion method of primal heuristic */
450  	   SCIP_CALL( SCIPincludeHeurOfins(scip) );
451  	
452  	   return SCIP_OKAY;
453  	}
454  	
455  	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
456  	static
457  	SCIP_DECL_HEURFREE(heurFreeOfins)
458  	{  /*lint --e{715}*/
459  	   SCIP_HEURDATA* heurdata;
460  	
461  	   assert(heur != NULL);
462  	   assert(scip != NULL);
463  	
464  	   /* get heuristic data */
465  	   heurdata = SCIPheurGetData(heur);
466  	   assert(heurdata != NULL);
467  	
468  	   /* free heuristic data */
469  	   SCIPfreeBlockMemory(scip, &heurdata);
470  	   SCIPheurSetData(heur, NULL);
471  	
472  	   return SCIP_OKAY;
473  	}
474  	
475  	/** execution method of primal heuristic */
476  	static
477  	SCIP_DECL_HEUREXEC(heurExecOfins)
478  	{/*lint --e{715}*/
479  	   SCIP_HEURDATA* heurdata;
480  	   SCIP_VAR** vars;
481  	   SCIP_Bool* chgcoeffs;
482  	   SCIP_Longint nstallnodes;
483  	   int nchgcoefs;
484  	   int nvars;
485  	   int v;
486  	
487  	   assert( heur != NULL );
488  	   assert( scip != NULL );
489  	   assert( result != NULL );
490  	
491  	   *result = SCIP_DELAYED;
492  	
493  	   /* do not call heuristic of node was already detected to be infeasible */
494  	   if( nodeinfeasible )
495  	      return SCIP_OKAY;
496  	
497  	   /* get heuristic data */
498  	   heurdata = SCIPheurGetData(heur);
499  	   assert( heurdata != NULL );
500  	
501  	   /* only call heuristic, if reoptimization is enabled */
502  	   if( !SCIPisReoptEnabled(scip) )
503  	      return SCIP_OKAY;
504  	
505  	   /* only call the heuristic, if we are in run >= 2 */
506  	   if( SCIPgetNReoptRuns(scip) <= 1 )
507  	      return SCIP_OKAY;
508  	
509  	   *result = SCIP_DIDNOTRUN;
510  	
511  	   if( SCIPisStopped(scip) )
512  	      return SCIP_OKAY;
513  	
514  	   /* calculate the maximal number of branching nodes until heuristic is aborted */
515  	   nstallnodes = (SCIP_Longint)(heurdata->nodesquot * SCIPgetNNodes(scip));
516  	
517  	   /* reward OFINS if it succeeded often */
518  	   nstallnodes = (SCIP_Longint)(nstallnodes * 3.0 * (SCIPheurGetNBestSolsFound(heur)+1.0)/(SCIPheurGetNCalls(heur) + 1.0));
519  	   nstallnodes -= 100 * SCIPheurGetNCalls(heur);  /* count the setup costs for the sub-SCIP as 100 nodes */
520  	   nstallnodes += heurdata->nodesofs;
521  	
522  	   /* determine the node limit for the current process */
523  	   nstallnodes = MIN(nstallnodes, heurdata->maxnodes);
524  	
525  	   /* check whether we have enough nodes left to call subproblem solving */
526  	   if( nstallnodes < heurdata->minnodes )
527  	   {
528  	      SCIPdebugMsg(scip, "skipping OFINS: nstallnodes=%" SCIP_LONGINT_FORMAT ", minnodes=%" SCIP_LONGINT_FORMAT "\n", nstallnodes, heurdata->minnodes);
529  	      return SCIP_OKAY;
530  	   }
531  	
532  	   /* get variable data and check which coefficient has changed  */
533  	   vars = SCIPgetVars(scip);
534  	   nvars = SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) + SCIPgetNImplVars(scip);
535  	   nchgcoefs = 0;
536  	
537  	   SCIP_CALL( SCIPallocBufferArray(scip, &chgcoeffs, nvars) );
538  	
539  	   for( v = 0; v < nvars; v++ )
540  	   {
541  	      SCIP_Real newcoef;
542  	      SCIP_Real oldcoef;
543  	      SCIP_Real newcoefabs;
544  	      SCIP_Real oldcoefabs;
545  	      SCIP_Real frac;
546  	
547  	      /* we only want to count variables that are unfixed after the presolving */
548  	      assert(SCIPvarGetStatus(vars[v]) != SCIP_VARSTATUS_ORIGINAL);
549  	      assert(SCIPvarIsActive(vars[v]));
550  	
551  	      SCIP_CALL( SCIPgetReoptOldObjCoef(scip, vars[v], SCIPgetNReoptRuns(scip), &newcoef) );
552  	      SCIP_CALL( SCIPgetReoptOldObjCoef(scip, vars[v], SCIPgetNReoptRuns(scip)-1, &oldcoef) );
553  	      newcoefabs = REALABS(newcoef);
554  	      oldcoefabs = REALABS(oldcoef);
555  	
556  	      /* if both coefficients are zero nothing has changed */
557  	      if( SCIPisZero(scip, newcoef) && SCIPisZero(scip, oldcoef) )
558  	      {
559  	         frac = 0;
560  	      }
561  	      /* if exactly one coefficient is zero, the other need to be close to zero */
562  	      else if( SCIPisZero(scip, newcoef) || SCIPisZero(scip, oldcoef) )
563  	      {
564  	         assert(SCIPisZero(scip, newcoef) != SCIPisZero(scip, oldcoef));
565  	         if( !SCIPisZero(scip, newcoef) )
566  	            frac = MIN(1, newcoefabs);
567  	         else
568  	            frac = MIN(1, oldcoefabs);
569  	      }
570  	      /* if both coefficients have the same sign we calculate the quotient
571  	       * MIN(newcoefabs, oldcoefabs)/MAX(newcoefabs, oldcoefabs)
572  	       */
573  	      else if( SCIPisPositive(scip, newcoef) == SCIPisPositive(scip, oldcoef) )
574  	      {
575  	         frac = 1.0 - MIN(newcoefabs, oldcoefabs)/MAX(newcoefabs, oldcoefabs);
576  	      }
577  	      /* if both coefficients have a different sign, we set frac = 1 */
578  	      else
579  	      {
580  	         assert((SCIPisPositive(scip, newcoef) && SCIPisNegative(scip, oldcoef))
581  	             || (SCIPisNegative(scip, newcoef) && SCIPisPositive(scip, oldcoef)));
582  	
583  	         frac = 1;
584  	      }
585  	
586  	      if( frac > heurdata->maxchange )
587  	      {
588  	         chgcoeffs[v] = TRUE;
589  	         nchgcoefs++;
590  	      }
591  	      else
592  	         chgcoeffs[v] = FALSE;
593  	   }
594  	
595  	   SCIPdebugMsg(scip, "%d (rate %.4f) changed coefficients\n", nchgcoefs, nchgcoefs/((SCIP_Real)nvars));
596  	
597  	   /* we only want to run the heuristic, if there at least 3 changed coefficients.
598  	    * if the number of changed coefficients is 2 the trivialnegation heuristic will construct an
599  	    * optimal solution without solving a MIP.
600  	    */
601  	   if( nchgcoefs < 3 )
602  	      goto TERMINATE;
603  	
604  	   /* run the heuristic, if not too many coefficients have changed */
605  	   if( nchgcoefs/((SCIP_Real)nvars) > heurdata->maxchangerate )
606  	      goto TERMINATE;
607  	
608  	   SCIP_CALL( applyOfins(scip, heur, heurdata, result, nstallnodes, chgcoeffs) );
609  	
610  	  TERMINATE:
611  	   SCIPfreeBufferArray(scip, &chgcoeffs);
612  	
613  	   return SCIP_OKAY;
614  	}
615  	
616  	
617  	/*
618  	 * primal heuristic specific interface methods
619  	 */
620  	
621  	/** creates the ofins primal heuristic and includes it in SCIP */
622  	SCIP_RETCODE SCIPincludeHeurOfins(
623  	   SCIP*                 scip                /**< SCIP data structure */
624  	   )
625  	{
626  	   SCIP_HEURDATA* heurdata;
627  	   SCIP_HEUR* heur;
628  	
629  	   /* create ofins primal heuristic data */
630  	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
631  	   assert(heurdata != NULL);
632  	
633  	   /* include primal heuristic */
634  	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
635  	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
636  	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecOfins, heurdata) );
637  	
638  	   assert(heur != NULL);
639  	
640  	   /* set non fundamental callbacks via setter functions */
641  	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyOfins) );
642  	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeOfins) );
643  	
644  	   /* add ofins primal heuristic parameters */
645  	
646  	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/maxnodes",
647  	         "maximum number of nodes to regard in the subproblem",
648  	         &heurdata->maxnodes, TRUE, DEFAULT_MAXNODES, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
649  	
650  	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/minnodes",
651  	         "minimum number of nodes required to start the subproblem",
652  	         &heurdata->minnodes, TRUE, DEFAULT_MINNODES, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
653  	
654  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/maxchangerate",
655  	         "maximal rate of changed coefficients",
656  	         &heurdata->maxchangerate, FALSE, DEFAULT_MAXCHGRATE, 0.0, 1.0, NULL, NULL) );
657  	
658  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/maxchange",
659  	         "maximal rate of change per coefficient to get fixed",
660  	         &heurdata->maxchange, FALSE, DEFAULT_MAXCHANGE, 0.0, 1.0, NULL, NULL) );
661  	
662  	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/copycuts",
663  	         "should all active cuts from cutpool be copied to constraints in subproblem?",
664  	         &heurdata->copycuts, TRUE, DEFAULT_COPYCUTS, NULL, NULL) );
665  	
666  	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/addallsols",
667  	         "should all subproblem solutions be added to the original SCIP?",
668  	         &heurdata->addallsols, TRUE, DEFAULT_ADDALLSOLS, NULL, NULL) );
669  	
670  	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/nodesofs",
671  	         "number of nodes added to the contingent of the total nodes",
672  	         &heurdata->nodesofs, FALSE, DEFAULT_NODESOFS, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
673  	
674  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/nodesquot",
675  	         "contingent of sub problem nodes in relation to the number of nodes of the original problem",
676  	         &heurdata->nodesquot, FALSE, DEFAULT_NODESQUOT, 0.0, 1.0, NULL, NULL) );
677  	
678  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/minimprove",
679  	         "factor by which RENS should at least improve the incumbent",
680  	         &heurdata->minimprove, TRUE, DEFAULT_MINIMPROVE, 0.0, 1.0, NULL, NULL) );
681  	
682  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/lplimfac",
683  	         "factor by which the limit on the number of LP depends on the node limit",
684  	         &heurdata->lplimfac, TRUE, DEFAULT_LPLIMFAC, 1.0, SCIP_REAL_MAX, NULL, NULL) );
685  	
686  	   return SCIP_OKAY;
687  	}
688