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_trustregion.c
26   	 * @ingroup DEFPLUGINS_HEUR
27   	 * @brief  Large neighborhood search heuristic for Benders' decomposition based on trust region methods
28   	 * @author Stephen J. Maher
29   	 *
30   	 * The Trust Region heuristic draws upon trust region methods for solving optimization problems, especially in the
31   	 * context of Benders' decomposition. This heuristic has been developed to improve the heuristic performance of the
32   	 * Benders' decomposition algorithm within SCIP.
33   	 *
34   	 * The Trust Region heuristic copies the original SCIP instance and adds a constraint to penalize changes from the
35   	 * incumbent solution. Consider a problem that includes a set of binary variables \f$\mathcal{B}\f$. Given a feasible
36   	 * solution \f$\hat{x}\f$ to the original problem, we define the set \f$\mathcal{B}^{+}\f$ as the index set for the
37   	 * binary variables that are 1 in the input solution and \f$\mathcal{B}^{-}\f$ as the index set for binary variables
38   	 * that are 0. The trust region constraint, which is added to the sub-SCIP, is given by
39   	 *
40   	 * \f[
41   	 *    \sum_{i \in \mathcal{B}^{+}}(1 - x_{i}) + \sum_{i \in \mathcal{B}^{-}}x_{i} \le \theta
42   	 * \f]
43   	 *
44   	 * The variable \f$\theta\f$ measure the distance, in terms of the binary variables, of candidate solutions to the input
45   	 * solution.
46   	 *
47   	 * In addition, an upper bounding constraint is explicitly added to enforce a minimum improvement from the heuristic,
48   	 * given by \f$f(x) \le f(\hat{x}) - \epsilon\f$. The parameter \f$\epsilon \ge 0\f$ denotes the minimum improvement
49   	 * that must be achieved by the heuristic.
50   	 *
51   	 * The objective function is then modified to \f$f(x) + M\theta\f$, where \f$M\f$ is a parameter for penalizing the
52   	 * distance of solutions from the input solution \f$\hat{x}\f$.
53   	 *
54   	 * If a new incumbent solution is found by this heuristic, then the Trust Region heuristic is immediately
55   	 * re-executed with this new incumbent solution.
56   	 */
57   	
58   	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
59   	
60   	#include "blockmemshell/memory.h"
61   	#include "scip/cons_linear.h"
62   	#include "scip/heuristics.h"
63   	#include "scip/heur_trustregion.h"
64   	#include "scip/pub_event.h"
65   	#include "scip/pub_heur.h"
66   	#include "scip/pub_message.h"
67   	#include "scip/pub_misc.h"
68   	#include "scip/pub_sol.h"
69   	#include "scip/pub_var.h"
70   	#include "scip/scip_branch.h"
71   	#include "scip/scip_cons.h"
72   	#include "scip/scip_copy.h"
73   	#include "scip/scip_event.h"
74   	#include "scip/scip_general.h"
75   	#include "scip/scip_heur.h"
76   	#include "scip/scip_mem.h"
77   	#include "scip/scip_message.h"
78   	#include "scip/scip_nodesel.h"
79   	#include "scip/scip_numerics.h"
80   	#include "scip/scip_param.h"
81   	#include "scip/scip_prob.h"
82   	#include "scip/scip_sol.h"
83   	#include "scip/scip_solve.h"
84   	#include "scip/scip_solvingstats.h"
85   	#include "scip/scip_var.h"
86   	#include <string.h>
87   	
88   	#define HEUR_NAME             "trustregion"
89   	#define HEUR_DESC             "LNS heuristic for Benders' decomposition based on trust region methods"
90   	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_LNS
91   	#define HEUR_PRIORITY         -1102010
92   	#define HEUR_FREQ             -1
93   	#define HEUR_FREQOFS          0
94   	#define HEUR_MAXDEPTH         -1
95   	#define HEUR_TIMING           SCIP_HEURTIMING_AFTERNODE
96   	#define HEUR_USESSUBSCIP      TRUE  /**< does the heuristic use a secondary SCIP instance? */
97   	
98   	#define DEFAULT_MINBINVARS    10        /**< the minimum number of binary variables necessary to run the heuristic */
99   	#define DEFAULT_NODESOFS      1000      /**< number of nodes added to the contingent of the total nodes */
100  	#define DEFAULT_MAXNODES      10000     /**< maximum number of nodes to regard in the subproblem */
101  	#define DEFAULT_MINNODES      100       /**< minimum number of nodes required to start the subproblem */
102  	#define DEFAULT_NODESQUOT     0.05      /**< contingent of sub problem nodes in relation to original nodes */
103  	#define DEFAULT_LPLIMFAC      1.5       /**< factor by which the limit on the number of LP depends on the node limit */
104  	#define DEFAULT_NWAITINGNODES 1         /**< number of nodes without incumbent change that heuristic should wait */
105  	#define DEFAULT_USELPROWS     FALSE     /**< should subproblem be created out of the rows in the LP rows,
106  	                                         * otherwise, the copy constructors of the constraints handlers are used */
107  	#define DEFAULT_COPYCUTS      TRUE      /**< if DEFAULT_USELPROWS is FALSE, then should all active cuts from the cutpool
108  	                                         * of the original scip be copied to constraints of the subscip */
109  	#define DEFAULT_BESTSOLLIMIT   3         /**< limit on number of improving incumbent solutions in sub-CIP */
110  	
111  	#define DEFAULT_VIOLPENALTY   100.0     /**< the penalty for violating the trust region */
112  	#define DEFAULT_OBJMINIMPROVE 1e-2      /**< the minimum absolute improvement in the objective function value */
113  	
114  	/* event handler properties */
115  	#define EVENTHDLR_NAME         "Trustregion"
116  	#define EVENTHDLR_DESC         "LP event handler for " HEUR_NAME " heuristic"
117  	
118  	
119  	#define EXECUTE               0
120  	#define WAITFORNEWSOL         1
121  	
122  	
123  	/*
124  	 * Data structures
125  	 */
126  	
127  	/** primal heuristic data */
128  	struct SCIP_HeurData
129  	{
130  	   SCIP_SOL*             lastsol;            /**< the last incumbent trustregion used as reference point */
131  	   SCIP_Longint          usednodes;          /**< amount of nodes trust region used during all calls */
132  	   SCIP_Real             nodesquot;          /**< contingent of sub problem nodes in relation to original nodes */
133  	   SCIP_Real             nodelimit;          /**< the nodelimit employed in the current sub-SCIP, for the event handler*/
134  	   SCIP_Real             lplimfac;           /**< factor by which the limit on the number of LP depends on the node limit */
135  	   SCIP_Real             violpenalty;        /**< the penalty for violating the trust region */
136  	   SCIP_Real             objminimprove;      /**< the minimum absolute improvement in the objective function value */
137  	   int                   nwaitingnodes;      /**< number of nodes without incumbent change that heuristic should wait */
138  	   int                   nodesofs;           /**< number of nodes added to the contingent of the total nodes */
139  	   int                   minnodes;           /**< minimum number of nodes required to start the subproblem */
140  	   int                   maxnodes;           /**< maximum number of nodes to regard in the subproblem */
141  	   int                   minbinvars;         /**< minimum number of binary variables necessary to run the heuristic */
142  	   int                   callstatus;         /**< current status of trustregion heuristic */
143  	   int                   curminnodes;        /**< current minimal number of nodes required to start the subproblem */
144  	   int                   bestsollimit;       /**< limit on number of improving incumbent solutions in sub-CIP */
145  	   SCIP_Bool             uselprows;          /**< should subproblem be created out of the rows in the LP rows? */
146  	   SCIP_Bool             copycuts;           /**< if uselprows == FALSE, should all active cuts from cutpool be copied
147  	                                              *   to constraints in subproblem? */
148  	};
149  	
150  	
151  	/*
152  	 * Local methods
153  	 */
154  	
155  	/** create the extra constraint of trust region and add it to \p subscip */
156  	static
157  	SCIP_RETCODE addTrustRegionConstraints(
158  	   SCIP*                 scip,               /**< SCIP data structure of the original problem */
159  	   SCIP*                 subscip,            /**< SCIP data structure of the subproblem */
160  	   SCIP_VAR**            subvars,            /**< variables of the subproblem */
161  	   SCIP_HEURDATA*        heurdata            /**< heuristic's data structure */
162  	   )
163  	{
164  	   SCIP_CONS* cons;                        /* trust region constraint to create */
165  	   SCIP_VAR** consvars;
166  	   SCIP_VAR** vars;
167  	   SCIP_SOL* bestsol;
168  	
169  	   int nvars;
170  	   int nbinvars;
171  	   int nconsvars;
172  	   int i;
173  	   SCIP_Real lhs;
174  	   SCIP_Real rhs;
175  	   SCIP_Real* consvals;
176  	   char name[SCIP_MAXSTRLEN];
177  	
178  	   /* adding the neighborhood constraint for the trust region heuristic */
179  	   SCIP_CALL( SCIPaddTrustregionNeighborhoodConstraint(scip, subscip, subvars, heurdata->violpenalty) );
180  	
181  	   /* get the data of the variables and the best solution */
182  	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, &nbinvars, NULL, NULL, NULL) );
183  	   bestsol = SCIPgetBestSol(scip);
184  	   assert( bestsol != NULL );
185  	
186  	   /* memory allocation */
187  	   SCIP_CALL( SCIPallocBufferArray(scip, &consvars, nvars + 1) );
188  	   SCIP_CALL( SCIPallocBufferArray(scip, &consvals, nvars + 1) );
189  	   nconsvars = 0;
190  	
191  	   /* create the upper bounding constraint. An absolute minimum improvement is used for this heuristic. This is
192  	    * different to other LNS heuristics, where a relative improvement is used. The absolute improvement tries to take
193  	    * into account problem specific information that is available to the user, such as a minimum step in the objective
194  	    * limit if the objective function is integer
195  	    */
196  	   lhs = -SCIPinfinity(subscip);
197  	   rhs = SCIPgetSolTransObj(scip, bestsol) - heurdata->objminimprove;
198  	
199  	   /* if the objective function is integer, then the floor of the RHS is taken */
200  	   if( SCIPisObjIntegral(scip) )
201  	      rhs = SCIPfeasFloor(scip, rhs);
202  	
203  	   /* adding the coefficients to the upper bounding constraint */
204  	   for( i = 0; i < nvars; i++ )
205  	   {
206  	      if( subvars[i] == NULL )
207  	         continue;
208  	      consvals[nconsvars] = SCIPvarGetObj(subvars[i]);
209  	      consvars[nconsvars] = subvars[i];
210  	      ++nconsvars;
211  	   }
212  	
213  	   /* creates trustregion constraint and adds it to subscip */
214  	   (void) SCIPsnprintf(name, SCIP_MAXSTRLEN, "%s_upperboundcons", SCIPgetProbName(scip));
215  	
216  	   SCIP_CALL( SCIPcreateConsLinear(subscip, &cons, name, nconsvars, consvars, consvals,
217  	         lhs, rhs, TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
218  	   SCIP_CALL( SCIPaddCons(subscip, cons) );
219  	   SCIP_CALL( SCIPreleaseCons(subscip, &cons) );
220  	
221  	   /* free local memory */
222  	   SCIPfreeBufferArray(scip, &consvals);
223  	   SCIPfreeBufferArray(scip, &consvars);
224  	
225  	   return SCIP_OKAY;
226  	}
227  	
228  	
229  	/* ---------------- Callback methods of event handler ---------------- */
230  	
231  	/** event handler execution callback to interrupt the solution process */
232  	static
233  	SCIP_DECL_EVENTEXEC(eventExecTrustregion)
234  	{
235  	   SCIP_HEURDATA* heurdata;
236  	
237  	   assert(eventhdlr != NULL);
238  	   assert(eventdata != NULL);
239  	   assert(strcmp(SCIPeventhdlrGetName(eventhdlr), EVENTHDLR_NAME) == 0);
240  	   assert(event != NULL);
241  	   assert(SCIPeventGetType(event) & SCIP_EVENTTYPE_LPSOLVED);
242  	
243  	   heurdata = (SCIP_HEURDATA*)eventdata;
244  	   assert(heurdata != NULL);
245  	
246  	   /* interrupt solution process of sub-SCIP */
247  	   if( SCIPgetNLPs(scip) > heurdata->lplimfac * heurdata->nodelimit )
248  	   {
249  	      SCIPdebugMsg(scip, "interrupt after  %" SCIP_LONGINT_FORMAT " LPs\n",SCIPgetNLPs(scip));
250  	      SCIP_CALL( SCIPinterruptSolve(scip) );
251  	   }
252  	
253  	   return SCIP_OKAY;
254  	}
255  	
256  	
257  	/*
258  	 * Callback methods of primal heuristic
259  	 */
260  	
261  	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
262  	static
263  	SCIP_DECL_HEURCOPY(heurCopyTrustregion)
264  	{  /*lint --e{715}*/
265  	   assert(scip != NULL);
266  	   assert(heur != NULL);
267  	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
268  	
269  	   /* call inclusion method of primal heuristic */
270  	   SCIP_CALL( SCIPincludeHeurTrustregion(scip) );
271  	
272  	   return SCIP_OKAY;
273  	}
274  	
275  	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
276  	static
277  	SCIP_DECL_HEURFREE(heurFreeTrustregion)
278  	{  /*lint --e{715}*/
279  	   SCIP_HEURDATA* heurdata;
280  	
281  	   assert( heur != NULL );
282  	   assert( scip != NULL );
283  	
284  	   /* get heuristic data */
285  	   heurdata = SCIPheurGetData(heur);
286  	   assert( heurdata != NULL );
287  	
288  	   /* free heuristic data */
289  	   SCIPfreeBlockMemory(scip, &heurdata);
290  	   SCIPheurSetData(heur, NULL);
291  	
292  	   return SCIP_OKAY;
293  	}
294  	
295  	
296  	/** initialization method of primal heuristic (called after problem was transformed) */
297  	static
298  	SCIP_DECL_HEURINIT(heurInitTrustregion)
299  	{  /*lint --e{715}*/
300  	   SCIP_HEURDATA* heurdata;
301  	
302  	   assert( heur != NULL );
303  	   assert( scip != NULL );
304  	
305  	   /* get heuristic's data */
306  	   heurdata = SCIPheurGetData(heur);
307  	   assert( heurdata != NULL );
308  	
309  	   /* with a little abuse we initialize the heurdata as if trustregion would have finished its last step regularly */
310  	   heurdata->callstatus = WAITFORNEWSOL;
311  	   heurdata->lastsol = NULL;
312  	   heurdata->usednodes = 0;
313  	   heurdata->curminnodes = heurdata->minnodes;
314  	
315  	   return SCIP_OKAY;
316  	}
317  	
318  	/** sets up and solves the sub SCIP for the Trust Region heuristic */
319  	static
320  	SCIP_RETCODE setupAndSolveSubscipTrustregion(
321  	   SCIP*                 scip,               /**< SCIP data structure */
322  	   SCIP*                 subscip,            /**< the subproblem created by trustregion */
323  	   SCIP_HEUR*            heur,               /**< trustregion heuristic */
324  	   SCIP_Longint          nsubnodes,          /**< nodelimit for subscip */
325  	   SCIP_RESULT*          result              /**< result pointer */
326  	   )
327  	{
328  	   SCIP_VAR** subvars;
329  	   SCIP_EVENTHDLR* eventhdlr;
330  	   SCIP_HEURDATA* heurdata;
331  	   SCIP_HASHMAP* varmapfw;
332  	   SCIP_VAR** vars;
333  	
334  	   int nvars;
335  	   int i;
336  	
337  	   SCIP_Bool success;
338  	
339  	   assert(scip != NULL);
340  	   assert(subscip != NULL);
341  	   assert(heur != NULL);
342  	
343  	   heurdata = SCIPheurGetData(heur);
344  	   assert(heurdata != NULL);
345  	
346  	   /* get the data of the variables and the best solution */
347  	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
348  	
349  	   /* create the variable mapping hash map */
350  	   SCIP_CALL( SCIPhashmapCreate(&varmapfw, SCIPblkmem(subscip), nvars) );
351  	   success = FALSE;
352  	
353  	   /* create a problem copy as sub SCIP */
354  	   SCIP_CALL( SCIPcopyLargeNeighborhoodSearch(scip, subscip, varmapfw, "trustregion", NULL, NULL, 0, heurdata->uselprows,
355  	         heurdata->copycuts, &success, NULL) );
356  	
357  	   SCIPdebugMsg(scip, "Copying SCIP was %s successful.\n", success ? "" : "not ");
358  	
359  	   /* if the subproblem could not be created, free memory and return */
360  	   if( !success )
361  	   {
362  	      *result = SCIP_DIDNOTRUN;
363  	      goto TERMINATE;
364  	   }
365  	
366  	   /* create event handler for LP events */
367  	   eventhdlr = NULL;
368  	   SCIP_CALL( SCIPincludeEventhdlrBasic(subscip, &eventhdlr, EVENTHDLR_NAME, EVENTHDLR_DESC, eventExecTrustregion, NULL) );
369  	   if( eventhdlr == NULL )
370  	   {
371  	      /* free hash map */
372  	      SCIPhashmapFree(&varmapfw);
373  	
374  	      SCIPerrorMessage("event handler for " HEUR_NAME " heuristic not found.\n");
375  	      return SCIP_PLUGINNOTFOUND;
376  	   }
377  	
378  	   SCIP_CALL( SCIPallocBufferArray(scip, &subvars, nvars) );
379  	   for (i = 0; i < nvars; ++i)
380  	      subvars[i] = (SCIP_VAR*) SCIPhashmapGetImage(varmapfw, vars[i]);
381  	
382  	   /* free hash map */
383  	   SCIPhashmapFree(&varmapfw);
384  	
385  	   heurdata->nodelimit = nsubnodes;
386  	   SCIP_CALL( SCIPsetCommonSubscipParams(scip, subscip, nsubnodes, MAX(10, nsubnodes/10), heurdata->bestsollimit) );
387  	
388  	   SCIP_CALL( addTrustRegionConstraints(scip, subscip, subvars, heurdata) );
389  	
390  	   /* catch LP events of sub-SCIP */
391  	   if( !heurdata->uselprows )
392  	   {
393  	      assert(eventhdlr != NULL);
394  	
395  	      SCIP_CALL( SCIPtransformProb(subscip) );
396  	      SCIP_CALL( SCIPcatchEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, NULL) );
397  	   }
398  	
399  	   /* solve the subproblem */
400  	   SCIPdebugMsg(scip, "solving trust region subproblem with maxnodes %" SCIP_LONGINT_FORMAT "\n", nsubnodes);
401  	
402  	   SCIP_CALL( SCIPsetIntParam(subscip, "heuristics/trysol/priority", 100000) );
403  	
404  	   /* Errors in solving the subproblem should not kill the overall solving process
405  	    * Hence, the return code is caught and a warning is printed, only in debug mode, SCIP will stop.
406  	    */
407  	   SCIP_CALL_ABORT( SCIPsolve(subscip) );
408  	
409  	   /* drop LP events of sub-SCIP */
410  	   if( !heurdata->uselprows )
411  	   {
412  	      assert(eventhdlr != NULL);
413  	
414  	      SCIP_CALL( SCIPdropEvent(subscip, SCIP_EVENTTYPE_LPSOLVED, eventhdlr, (SCIP_EVENTDATA*) heurdata, -1) );
415  	   }
416  	
417  	   /* print solving statistics of subproblem if we are in SCIP's debug mode */
418  	   SCIPdebug( SCIP_CALL( SCIPprintStatistics(subscip, NULL) ) );
419  	
420  	   heurdata->usednodes += SCIPgetNNodes(subscip);
421  	   SCIPdebugMsg(scip, "trust region used %" SCIP_LONGINT_FORMAT "/%" SCIP_LONGINT_FORMAT " nodes\n",
422  	      SCIPgetNNodes(subscip), nsubnodes);
423  	
424  	   /* checks the solutions of the sub SCIP and adds them to the main SCIP if feasible */
425  	   SCIP_CALL( SCIPtranslateSubSols(scip, subscip, heur, subvars, &success, NULL) );
426  	
427  	   if( success )
428  	      *result = SCIP_FOUNDSOL;
429  	
430  	   /* checking the status of the subscip */
431  	   heurdata->callstatus = WAITFORNEWSOL;
432  	   if( SCIPgetStatus(subscip) == SCIP_STATUS_NODELIMIT || SCIPgetStatus(subscip) == SCIP_STATUS_STALLNODELIMIT
433  	      || SCIPgetStatus(subscip) == SCIP_STATUS_TOTALNODELIMIT )
434  	   {
435  	      heurdata->callstatus = EXECUTE;
436  	      heurdata->curminnodes *= 2;
437  	   }
438  	
439  	 TERMINATE:
440  	   /* free subproblem */
441  	   SCIPfreeBufferArray(scip, &subvars);
442  	
443  	   return SCIP_OKAY;
444  	}
445  	
446  	
447  	/** execution method of primal heuristic */
448  	static
449  	SCIP_DECL_HEUREXEC(heurExecTrustregion)
450  	{  /*lint --e{715}*/
451  	   SCIP_Longint maxnnodes;
452  	   SCIP_Longint nsubnodes;
453  	
454  	   SCIP_HEURDATA* heurdata;
455  	   SCIP* subscip;
456  	
457  	   SCIP_SOL* bestsol;
458  	
459  	   SCIP_Bool success;
460  	   SCIP_RETCODE retcode;
461  	
462  	   assert(heur != NULL);
463  	   assert(scip != NULL);
464  	   assert(result != NULL);
465  	
466  	   *result = SCIP_DIDNOTRUN;
467  	
468  	   /* get heuristic's data */
469  	   heurdata = SCIPheurGetData(heur);
470  	   assert( heurdata != NULL );
471  	
472  	   /* there should be enough binary variables that a trust region constraint makes sense */
473  	   if( SCIPgetNBinVars(scip) < heurdata->minbinvars )
474  	      return SCIP_OKAY;
475  	
476  	   *result = SCIP_DELAYED;
477  	
478  	   /* only call heuristic, if an IP solution is at hand */
479  	   if( SCIPgetNSols(scip) <= 0  )
480  	      return SCIP_OKAY;
481  	
482  	   bestsol = SCIPgetBestSol(scip);
483  	   assert(bestsol != NULL);
484  	
485  	   /* only call heuristic, if the best solution comes from transformed problem */
486  	   if( SCIPsolIsOriginal(bestsol) )
487  	      return SCIP_OKAY;
488  	
489  	   /* only call heuristic, if enough nodes were processed since last incumbent */
490  	   if( SCIPgetNNodes(scip) - SCIPgetSolNodenum(scip, bestsol)  < heurdata->nwaitingnodes)
491  	      return SCIP_OKAY;
492  	
493  	   /* only call heuristic, if the best solution does not come from trivial heuristic */
494  	   if( SCIPsolGetHeur(bestsol) != NULL && strcmp(SCIPheurGetName(SCIPsolGetHeur(bestsol)), "trivial") == 0 )
495  	      return SCIP_OKAY;
496  	
497  	   /* calculate the maximal number of branching nodes until heuristic is aborted */
498  	   maxnnodes = (SCIP_Longint)(heurdata->nodesquot * SCIPgetNNodes(scip));
499  	
500  	   /* reward trust region if it found solutions often.
501  	    * In this case, the trust region heuristic is designed for Benders' decomposition and solutions found may not be
502  	    * added by this heuristic but by trysol. So we don't reward finding best solutions, but finding any solution. */
503  	   maxnnodes = (SCIP_Longint)(maxnnodes * (1.0 + 2.0*(SCIPheurGetNSolsFound(heur)+1.0)/(SCIPheurGetNCalls(heur)+1.0)));
504  	   maxnnodes -= 100 * SCIPheurGetNCalls(heur);  /* count the setup costs for the sub-MIP as 100 nodes */
505  	   maxnnodes += heurdata->nodesofs;
506  	
507  	   *result = SCIP_DIDNOTRUN;
508  	
509  	   /* we continue to execute the trust region heuristic until no new best solution is found */
510  	   do
511  	   {
512  	      SCIP_RESULT heurresult;
513  	
514  	      /* storing the best solution again since it is needed for the execution loop */
515  	      bestsol = SCIPgetBestSol(scip);
516  	
517  	      /* reset minnodes if new solution was found */
518  	      if( heurdata->lastsol != bestsol )
519  	      {
520  	         heurdata->curminnodes = heurdata->minnodes;
521  	         heurdata->callstatus = EXECUTE;
522  	         heurdata->lastsol = bestsol;
523  	      }
524  	
525  	      /* if no new solution was found and trust region also seems to fail, just keep on waiting */
526  	      if( heurdata->callstatus == WAITFORNEWSOL )
527  	         return SCIP_OKAY;
528  	
529  	      /* determine the node limit for the current process */
530  	      nsubnodes = maxnnodes - heurdata->usednodes;
531  	      nsubnodes = MIN(nsubnodes, heurdata->maxnodes);
532  	
533  	      /* check whether we have enough nodes left to call sub problem solving */
534  	      if( nsubnodes < heurdata->curminnodes )
535  	         return SCIP_OKAY;
536  	
537  	      if( SCIPisStopped(scip) )
538  	         return SCIP_OKAY;
539  	
540  	      /* check whether there is enough time and memory left */
541  	      SCIP_CALL( SCIPcheckCopyLimits(scip, &success) );
542  	
543  	      /* abort if no time is left or there is not enough memory to create a copy of SCIP */
544  	      if( !success )
545  	         return SCIP_OKAY;
546  	
547  	      heurresult = SCIP_DIDNOTFIND;
548  	
549  	      SCIPdebugMsg(scip, "running trust region heuristic ...\n");
550  	
551  	      SCIP_CALL( SCIPcreate(&subscip) );
552  	
553  	      retcode = setupAndSolveSubscipTrustregion(scip, subscip, heur, nsubnodes, &heurresult);
554  	
555  	      SCIP_CALL( SCIPfree(&subscip) );
556  	
557  	      /* if the result is FOUNDSOL, this means that a solution was found during a previous execution of the heuristic.
558  	       * So the heuristic result should only be updated if the result is not FOUNDSOL.
559  	       */
560  	      if( *result != SCIP_FOUNDSOL )
561  	         *result = heurresult;
562  	   }
563  	   while( bestsol != SCIPgetBestSol(scip) && retcode == SCIP_OKAY );
564  	
565  	   return retcode;
566  	}
567  	
568  	
569  	/*
570  	 * primal heuristic specific interface methods
571  	 */
572  	
573  	/** creates the trustregion primal heuristic and includes it in SCIP */
574  	SCIP_RETCODE SCIPincludeHeurTrustregion(
575  	   SCIP*                 scip                /**< SCIP data structure */
576  	   )
577  	{
578  	   SCIP_HEURDATA* heurdata;
579  	   SCIP_HEUR* heur;
580  	
581  	   /* create Trustregion primal heuristic data */
582  	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
583  	
584  	   /* include primal heuristic */
585  	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
586  	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
587  	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecTrustregion, heurdata) );
588  	
589  	   assert(heur != NULL);
590  	
591  	   /* set non-NULL pointers to callback methods */
592  	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyTrustregion) );
593  	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeTrustregion) );
594  	   SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitTrustregion) );
595  	
596  	   /* add trustregion primal heuristic parameters */
597  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/nodesofs",
598  	         "number of nodes added to the contingent of the total nodes",
599  	         &heurdata->nodesofs, FALSE, DEFAULT_NODESOFS, 0, INT_MAX, NULL, NULL) );
600  	
601  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/minbinvars",
602  	         "the number of binary variables necessary to run the heuristic",
603  	         &heurdata->minbinvars, FALSE, DEFAULT_MINBINVARS, 1, INT_MAX, NULL, NULL) );
604  	
605  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/nodesquot",
606  	         "contingent of sub problem nodes in relation to the number of nodes of the original problem",
607  	         &heurdata->nodesquot, FALSE, DEFAULT_NODESQUOT, 0.0, 1.0, NULL, NULL) );
608  	
609  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/lplimfac",
610  	         "factor by which the limit on the number of LP depends on the node limit",
611  	         &heurdata->lplimfac, TRUE, DEFAULT_LPLIMFAC, 1.0, SCIP_REAL_MAX, NULL, NULL) );
612  	
613  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/minnodes",
614  	         "minimum number of nodes required to start the subproblem",
615  	         &heurdata->minnodes, TRUE, DEFAULT_MINNODES, 0, INT_MAX, NULL, NULL) );
616  	
617  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/maxnodes",
618  	         "maximum number of nodes to regard in the subproblem",
619  	         &heurdata->maxnodes, TRUE, DEFAULT_MAXNODES, 0, INT_MAX, NULL, NULL) );
620  	
621  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/nwaitingnodes",
622  	         "number of nodes without incumbent change that heuristic should wait",
623  	         &heurdata->nwaitingnodes, TRUE, DEFAULT_NWAITINGNODES, 0, INT_MAX, NULL, NULL) );
624  	
625  	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/uselprows",
626  	         "should subproblem be created out of the rows in the LP rows?",
627  	         &heurdata->uselprows, TRUE, DEFAULT_USELPROWS, NULL, NULL) );
628  	
629  	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/copycuts",
630  	         "if uselprows == FALSE, should all active cuts from cutpool be copied to constraints in subproblem?",
631  	         &heurdata->copycuts, TRUE, DEFAULT_COPYCUTS, NULL, NULL) );
632  	
633  	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/bestsollimit",
634  	         "limit on number of improving incumbent solutions in sub-CIP",
635  	         &heurdata->bestsollimit, FALSE, DEFAULT_BESTSOLLIMIT, -1, INT_MAX, NULL, NULL) );
636  	
637  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/violpenalty",
638  	         "the penalty for each change in the binary variables from the candidate solution",
639  	         &heurdata->violpenalty, FALSE, DEFAULT_VIOLPENALTY, 0.0, SCIP_REAL_MAX, NULL, NULL) );
640  	
641  	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/objminimprove",
642  	         "the minimum absolute improvement in the objective function value",
643  	         &heurdata->objminimprove, FALSE, DEFAULT_OBJMINIMPROVE, 0.0, SCIP_REAL_MAX, NULL, NULL) );
644  	
645  	   return SCIP_OKAY;
646  	}
647