/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
	/*                                                                           */
	/*                  This file is part of the program and library             */
	/*         SCIP --- Solving Constraint Integer Programs                      */
	/*                                                                           */
	/*    Copyright (C) 2002-2021 Konrad-Zuse-Zentrum                            */
	/*                            fuer Informationstechnik Berlin                */
	/*                                                                           */
	/*  SCIP is distributed under the terms of the ZIB Academic License.         */
	/*                                                                           */
	/*  You should have received a copy of the ZIB Academic License              */
	/*  along with SCIP; see the file COPYING. If not visit scipopt.org.         */
	/*                                                                           */
	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
	
	/**@file   heur_alns.c
	 * @ingroup DEFPLUGINS_HEUR
	 * @brief  Adaptive large neighborhood search heuristic that orchestrates popular LNS heuristics
	 * @author Gregor Hendel
	 */
	
	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
	
	#include "blockmemshell/memory.h"
	#include "scip/cons_linear.h"
	#include "scip/heur_alns.h"
	#include "scip/heuristics.h"
	#include "scip/pub_bandit_epsgreedy.h"
	#include "scip/pub_bandit_exp3.h"
	#include "scip/pub_bandit.h"
	#include "scip/pub_bandit_ucb.h"
	#include "scip/pub_cons.h"
	#include "scip/pub_event.h"
	#include "scip/pub_heur.h"
	#include "scip/pub_message.h"
	#include "scip/pub_misc.h"
	#include "scip/pub_misc_select.h"
	#include "scip/pub_sol.h"
	#include "scip/pub_var.h"
	#include "scip/scip_bandit.h"
	#include "scip/scip_branch.h"
	#include "scip/scip_cons.h"
	#include "scip/scip_copy.h"
	#include "scip/scip_event.h"
	#include "scip/scip_general.h"
	#include "scip/scip_heur.h"
	#include "scip/scip_lp.h"
	#include "scip/scip_mem.h"
	#include "scip/scip_message.h"
	#include "scip/scip_nodesel.h"
	#include "scip/scip_numerics.h"
	#include "scip/scip_param.h"
	#include "scip/scip_prob.h"
	#include "scip/scip_randnumgen.h"
	#include "scip/scip_sol.h"
	#include "scip/scip_solve.h"
	#include "scip/scip_solvingstats.h"
	#include "scip/scip_table.h"
	#include "scip/scip_timing.h"
	#include "scip/scip_tree.h"
	#include "scip/scip_var.h"
	#include <string.h>
	
	#if !defined(_WIN32) && !defined(_WIN64)
	#include <strings.h> /*lint --e{766}*/ /* needed for strncasecmp() */
	#endif
	
	#define HEUR_NAME             "alns"
	#define HEUR_DESC             "Large neighborhood search heuristic that orchestrates the popular neighborhoods Local Branching, RINS, RENS, DINS etc."
	#define HEUR_DISPCHAR         SCIP_HEURDISPCHAR_LNS
	#define HEUR_PRIORITY         -1100500
	#define HEUR_FREQ             20
	#define HEUR_FREQOFS          0
	#define HEUR_MAXDEPTH         -1
	#define HEUR_TIMING           SCIP_HEURTIMING_AFTERNODE | SCIP_HEURTIMING_DURINGLPLOOP
	#define HEUR_USESSUBSCIP      TRUE  /**< does the heuristic use a secondary SCIP instance? */
	
	#define NNEIGHBORHOODS 9
	
	/*
	 * limit parameters for sub-SCIPs
	 */
	#define DEFAULT_NODESQUOT        0.1
	#define DEFAULT_NODESQUOTMIN     0.0
	#define DEFAULT_NODESOFFSET      500LL
	#define DEFAULT_NSOLSLIM         3
	#define DEFAULT_MINNODES         50LL
	#define DEFAULT_MAXNODES         5000LL
	#define DEFAULT_WAITINGNODES     25LL  /**< number of nodes since last incumbent solution that the heuristic should wait */
	#define DEFAULT_TARGETNODEFACTOR 1.05
	#define LRATEMIN                 0.01 /**<  lower bound for learning rate for target nodes and minimum improvement */
	#define LPLIMFAC                 4.0
	#define DEFAULT_INITDURINGROOT FALSE
	#define DEFAULT_MAXCALLSSAMESOL  -1   /**< number of allowed executions of the heuristic on the same incumbent solution */
	
	/*
	 * parameters for the minimum improvement
	 */
	#define DEFAULT_MINIMPROVELOW    0.01
	#define DEFAULT_MINIMPROVEHIGH   0.01
	#define MINIMPROVEFAC            1.5
	#define DEFAULT_STARTMINIMPROVE  0.01
	#define DEFAULT_ADJUSTMINIMPROVE FALSE
	#define DEFAULT_ADJUSTTARGETNODES TRUE /**< should the target nodes be dynamically adjusted? */
	
	/*
	 * bandit algorithm parameters
	 */
	#define DEFAULT_BESTSOLWEIGHT  1
	#define DEFAULT_BANDITALGO     'u'  /**< the default bandit algorithm: (u)pper confidence bounds, (e)xp.3, epsilon (g)reedy */
	#define DEFAULT_REWARDCONTROL  0.8  /**< reward control to increase the weight of the simple solution indicator and decrease the weight of the closed gap reward */
	#define DEFAULT_SCALEBYEFFORT  TRUE /**< should the reward be scaled by the effort? */
	#define DEFAULT_RESETWEIGHTS   TRUE /**< should the bandit algorithms be reset when a new problem is read? */
	#define DEFAULT_SUBSCIPRANDSEEDS FALSE /**< should random seeds of sub-SCIPs be altered to increase diversification? */
	#define DEFAULT_REWARDBASELINE 0.5  /**< the reward baseline to separate successful and failed calls */
	#define DEFAULT_FIXTOL         0.1  /**< tolerance by which the fixing rate may be missed without generic fixing */
	#define DEFAULT_UNFIXTOL       0.1  /**< tolerance by which the fixing rate may be exceeded without generic unfixing */
	#define DEFAULT_USELOCALREDCOST FALSE /**< should local reduced costs be used for generic (un)fixing? */
	#define DEFAULT_BETA           0.0  /**< default reward offset between 0 and 1 at every observation for exp3 */
	
	/*
	 * the following 3 parameters have been tuned by a simulation experiment
	 * as described in the paper.
	 */
	#define DEFAULT_EPS            0.4685844   /**< increase exploration in epsilon-greedy bandit algorithm */
	#define DEFAULT_ALPHA          0.0016      /**< parameter to increase the confidence width in UCB */
	#define DEFAULT_GAMMA          0.07041455  /**< default weight between uniform (gamma ~ 1) and weight driven (gamma ~ 0) probability distribution for exp3 */
	/*
	 * parameters to control variable fixing
	 */
	#define DEFAULT_USEREDCOST       TRUE  /**< should reduced cost scores be used for variable priorization? */
	#define DEFAULT_USEPSCOST        TRUE  /**< should pseudo cost scores be used for variable priorization? */
	#define DEFAULT_USEDISTANCES     TRUE  /**< should distances from fixed variables be used for variable priorization */
	#define DEFAULT_DOMOREFIXINGS    TRUE  /**< should the ALNS heuristic do more fixings by itself based on variable prioritization
	                                         *  until the target fixing rate is reached? */
	#define DEFAULT_ADJUSTFIXINGRATE TRUE  /**< should the heuristic adjust the target fixing rate based on the success? */
	#define FIXINGRATE_DECAY         0.75  /**< geometric decay for fixing rate adjustments */
	#define FIXINGRATE_STARTINC      0.2   /**< initial increment value for fixing rate */
	#define DEFAULT_USESUBSCIPHEURS  FALSE /**< should the heuristic activate other sub-SCIP heuristics during its search?  */
	#define DEFAULT_COPYCUTS         FALSE /**< should cutting planes be copied to the sub-SCIP? */
	#define DEFAULT_REWARDFILENAME   "-"   /**< file name to store all rewards and the selection of the bandit */
	
	/* individual random seeds */
	#define DEFAULT_SEED 113
	#define MUTATIONSEED 121
	#define CROSSOVERSEED 321
	
	/* individual neighborhood parameters */
	#define DEFAULT_MINFIXINGRATE_RENS 0.3
	#define DEFAULT_MAXFIXINGRATE_RENS 0.9
	#define DEFAULT_ACTIVE_RENS TRUE
	#define DEFAULT_PRIORITY_RENS 1.0
	
	#define DEFAULT_MINFIXINGRATE_RINS 0.3
	#define DEFAULT_MAXFIXINGRATE_RINS 0.9
	#define DEFAULT_ACTIVE_RINS TRUE
	#define DEFAULT_PRIORITY_RINS 1.0
	
	#define DEFAULT_MINFIXINGRATE_MUTATION 0.3
	#define DEFAULT_MAXFIXINGRATE_MUTATION 0.9
	#define DEFAULT_ACTIVE_MUTATION TRUE
	#define DEFAULT_PRIORITY_MUTATION 1.0
	
	#define DEFAULT_MINFIXINGRATE_LOCALBRANCHING 0.3
	#define DEFAULT_MAXFIXINGRATE_LOCALBRANCHING 0.9
	#define DEFAULT_ACTIVE_LOCALBRANCHING TRUE
	#define DEFAULT_PRIORITY_LOCALBRANCHING 1.0
	
	#define DEFAULT_MINFIXINGRATE_PROXIMITY 0.3
	#define DEFAULT_MAXFIXINGRATE_PROXIMITY 0.9
	#define DEFAULT_ACTIVE_PROXIMITY TRUE
	#define DEFAULT_PRIORITY_PROXIMITY 1.0
	
	#define DEFAULT_MINFIXINGRATE_CROSSOVER 0.3
	#define DEFAULT_MAXFIXINGRATE_CROSSOVER 0.9
	#define DEFAULT_ACTIVE_CROSSOVER TRUE
	#define DEFAULT_PRIORITY_CROSSOVER 1.0
	
	#define DEFAULT_MINFIXINGRATE_ZEROOBJECTIVE 0.3
	#define DEFAULT_MAXFIXINGRATE_ZEROOBJECTIVE 0.9
	#define DEFAULT_ACTIVE_ZEROOBJECTIVE TRUE
	#define DEFAULT_PRIORITY_ZEROOBJECTIVE 1.0
	
	#define DEFAULT_MINFIXINGRATE_DINS 0.3
	#define DEFAULT_MAXFIXINGRATE_DINS 0.9
	#define DEFAULT_ACTIVE_DINS TRUE
	#define DEFAULT_PRIORITY_DINS 1.0
	
	#define DEFAULT_MINFIXINGRATE_TRUSTREGION 0.3
	#define DEFAULT_MAXFIXINGRATE_TRUSTREGION 0.9
	#define DEFAULT_ACTIVE_TRUSTREGION FALSE
	#define DEFAULT_PRIORITY_TRUSTREGION 1.0
	
	
	#define DEFAULT_NSOLS_CROSSOVER 2 /**< parameter for the number of solutions that crossover should combine */
	#define DEFAULT_NPOOLSOLS_DINS  5 /**< number of pool solutions where binary solution values must agree */
	#define DEFAULT_VIOLPENALTY_TRUSTREGION 100.0  /**< the penalty for violating the trust region */
	
	/* event handler properties */
	#define EVENTHDLR_NAME         "Alns"
	#define EVENTHDLR_DESC         "LP event handler for " HEUR_NAME " heuristic"
	#define SCIP_EVENTTYPE_ALNS (SCIP_EVENTTYPE_LPSOLVED | SCIP_EVENTTYPE_SOLFOUND | SCIP_EVENTTYPE_BESTSOLFOUND)
	
	/* properties of the ALNS neighborhood statistics table */
	#define TABLE_NAME_NEIGHBORHOOD                  "neighborhood"
	#define TABLE_DESC_NEIGHBORHOOD                  "ALNS neighborhood statistics"
	#define TABLE_POSITION_NEIGHBORHOOD              12500                  /**< the position of the statistics table */
	#define TABLE_EARLIEST_STAGE_NEIGHBORHOOD        SCIP_STAGE_TRANSFORMED /**< output of the statistics table is only printed from this stage onwards */
	
	
	/** reward types of ALNS */
	enum RewardType {
	   REWARDTYPE_TOTAL,                         /**< combination of the other rewards */
	   REWARDTYPE_BESTSOL,                       /**< 1, if a new solution was found, 0 otherwise */
	   REWARDTYPE_CLOSEDGAP,                           /**< 0 if no solution was found, closed gap otherwise */
	   REWARDTYPE_NOSOLPENALTY,                  /**< 1 if a solution was found, otherwise between 0 and 1 depending on the effort spent  */
	   NREWARDTYPES
	};
	
	/*
	 * Data structures
	 */
	
	/*
	 * additional neighborhood data structures
	 */
	
	
	typedef struct data_crossover DATA_CROSSOVER; /**< crossover neighborhood data structure */
	
	typedef struct data_mutation DATA_MUTATION; /**< mutation neighborhood data structure */
	
	typedef struct data_dins DATA_DINS; /**< dins neighborhood data structure */
	
	typedef struct data_trustregion DATA_TRUSTREGION; /**< trustregion neighborhood data structure */
	
	typedef struct NH_FixingRate NH_FIXINGRATE; /** fixing rate data structure */
	
	typedef struct NH_Stats NH_STATS; /**< neighborhood statistics data structure */
	
	typedef struct Nh NH;             /**< neighborhood data structure */
	
	
	/*
	 * variable priorization data structure for sorting
	 */
	typedef struct VarPrio VARPRIO;
	
	/** callback to collect variable fixings of neighborhood */
	 #define DECL_VARFIXINGS(x) SCIP_RETCODE x ( \
	   SCIP*                 scip,               /**< SCIP data structure */                     \
	   NH*                   neighborhood,       /**< ALNS neighborhood data structure */         \
	   SCIP_VAR**            varbuf,             /**< buffer array to collect variables to fix */\
	   SCIP_Real*            valbuf,             /**< buffer array to collect fixing values */   \
	   int*                  nfixings,           /**< pointer to store the number of fixings */  \
	   SCIP_RESULT*          result              /**< result pointer */                          \
	   )
	
	/** callback for subproblem changes other than variable fixings
	 *
	 *  this callback can be used to further modify the subproblem by changes other than variable fixings.
	 *  Typical modifications include restrictions of variable domains, the formulation of additional constraints,
	 *  or changed objective coefficients.
	 *
	 *  The callback should set the \p success pointer to indicate whether it was successful with its modifications or not.
	 */
	#define DECL_CHANGESUBSCIP(x) SCIP_RETCODE x (  \
	   SCIP*                 sourcescip,         /**< source SCIP data structure */\
	   SCIP*                 targetscip,         /**< target SCIP data structure */\
	   NH*                   neighborhood,       /**< ALNS neighborhood data structure */\
	   SCIP_VAR**            subvars,            /**< array of targetscip variables in the same order as the source SCIP variables */\
	   int*                  ndomchgs,           /**< pointer to store the number of performed domain changes */\
	   int*                  nchgobjs,           /**< pointer to store the number of changed objective coefficients */ \
	   int*                  naddedconss,        /**< pointer to store the number of additional constraints */\
	   SCIP_Bool*            success             /**< pointer to store if the sub-MIP was successfully adjusted */\
	   )
	
	/** optional initialization callback for neighborhoods when a new problem is read */
	#define DECL_NHINIT(x) SCIP_RETCODE x (                                          \
	   SCIP*                 scip,               /**< SCIP data structure */         \
	   NH*                   neighborhood        /**< neighborhood data structure */ \
	   )
	
	/** deinitialization callback for neighborhoods when exiting a problem */
	#define DECL_NHEXIT(x) SCIP_RETCODE x ( \
	   SCIP*                 scip,               /**< SCIP data structure */         \
	   NH*                   neighborhood        /**< neighborhood data structure */ \
	   )
	
	/** deinitialization callback for neighborhoods before SCIP is freed */
	#define DECL_NHFREE(x) SCIP_RETCODE x (      \
	   SCIP*                 scip,               /**< SCIP data structure */         \
	   NH*                   neighborhood        /**< neighborhood data structure */ \
	   )
	
	/** callback function to return a feasible reference solution for further fixings
	 *
	 *  The reference solution should be stored in the \p solptr.
	 *  The \p result pointer can be used to indicate either
	 *
	 *  - SCIP_SUCCESS or
	 *  - SCIP_DIDNOTFIND
	 */
	#define DECL_NHREFSOL(x) SCIP_RETCODE x (                                       \
	   SCIP*                 scip,               /**< SCIP data structure */  \
	   NH*                   neighborhood,       /**< neighborhood data structure */ \
	   SCIP_SOL**            solptr,             /**< pointer to store the reference solution */ \
	   SCIP_RESULT*          result              /**< pointer to indicate the callback success whether a reference solution is available */ \
	   )
	
	/** callback function to deactivate neighborhoods on problems where they are irrelevant */
	#define DECL_NHDEACTIVATE(x) SCIP_RETCODE x (\
	   SCIP*                 scip,               /**< SCIP data structure */  \
	   SCIP_Bool*            deactivate          /**< pointer to store whether the neighborhood should be deactivated (TRUE) for an instance */ \
	   )
	
	/** sub-SCIP status code enumerator */
	enum HistIndex
	{
	   HIDX_OPT              = 0,                /**< sub-SCIP was solved to optimality  */
	   HIDX_USR              = 1,                /**< sub-SCIP was user interrupted */
	   HIDX_NODELIM          = 2,                /**< sub-SCIP reached the node limit */
	   HIDX_STALLNODE        = 3,                /**< sub-SCIP reached the stall node limit */
	   HIDX_INFEAS           = 4,                /**< sub-SCIP was infeasible */
	   HIDX_SOLLIM           = 5,                /**< sub-SCIP reached the solution limit */
	   HIDX_OTHER            = 6                 /**< sub-SCIP reached none of the above codes */
	};
	typedef enum HistIndex HISTINDEX;
	#define NHISTENTRIES 7
	
	
	/** statistics for a neighborhood */
	struct NH_Stats
	{
	   SCIP_CLOCK*           setupclock;         /**< clock for sub-SCIP setup time */
	   SCIP_CLOCK*           submipclock;        /**< clock for the sub-SCIP solve */
	   SCIP_Longint          usednodes;          /**< total number of used nodes */
	   SCIP_Real             oldupperbound;      /**< upper bound before the sub-SCIP started */
	   SCIP_Real             newupperbound;      /**< new upper bound for allrewards mode to work correctly */
	   int                   nruns;              /**< number of runs of a neighborhood */
	   int                   nrunsbestsol;       /**< number of runs that produced a new incumbent */
	   SCIP_Longint          nsolsfound;         /**< the total number of solutions found */
	   SCIP_Longint          nbestsolsfound;     /**< the total number of improving solutions found */
	   int                   nfixings;           /**< the number of fixings in one run */
	   int                   statushist[NHISTENTRIES]; /**< array to count sub-SCIP statuses */
	};
	
	
	/** fixing rate data structure to control the amount of target fixings of a neighborhood */
	struct NH_FixingRate
	{
	   SCIP_Real             minfixingrate;      /**< the minimum fixing rate */
	   SCIP_Real             targetfixingrate;   /**< the current target fixing rate */
	   SCIP_Real             increment;          /**< the current increment by which the target fixing rate is in-/decreased */
	   SCIP_Real             maxfixingrate;      /**< the maximum fixing rate */
	};
	
	/** neighborhood data structure with callbacks, statistics, fixing rate */
	struct Nh
	{
	   char*                 name;               /**< the name of this neighborhood */
	   NH_FIXINGRATE         fixingrate;         /**< fixing rate for this neighborhood */
	   NH_STATS              stats;              /**< statistics for this neighborhood */
	   DECL_VARFIXINGS       ((*varfixings));    /**< variable fixings callback for this neighborhood */
	   DECL_CHANGESUBSCIP    ((*changesubscip)); /**< callback for subproblem changes other than variable fixings */
	   DECL_NHINIT           ((*nhinit));        /**< initialization callback when a new problem is read */
	   DECL_NHEXIT           ((*nhexit));        /**< deinitialization callback when exiting a problem */
	   DECL_NHFREE           ((*nhfree));        /**< deinitialization callback before SCIP is freed */
	   DECL_NHREFSOL         ((*nhrefsol));      /**< callback function to return a reference solution for further fixings, or NULL */
	   DECL_NHDEACTIVATE     ((*nhdeactivate));  /**< callback function to deactivate neighborhoods on problems where they are irrelevant, or NULL if it is always active */
	   SCIP_Bool             active;             /**< is this neighborhood active or not? */
	   SCIP_Real             priority;           /**< positive call priority to initialize bandit algorithms */
	   union
	   {
	      DATA_MUTATION*     mutation;           /**< mutation data */
	      DATA_CROSSOVER*    crossover;          /**< crossover data */
	      DATA_DINS*         dins;               /**< dins data */
	      DATA_TRUSTREGION*  trustregion;        /**< trustregion data */
	   }                     data;               /**< data object for neighborhood specific data */
	};
	
	/** mutation neighborhood data structure */
	struct data_mutation
	{
	   SCIP_RANDNUMGEN*      rng;                /**< random number generator */
	};
	
	/** crossover neighborhood data structure */
	struct data_crossover
	{
	   int                   nsols;              /**< the number of solutions that crossover should combine */
	   SCIP_RANDNUMGEN*      rng;                /**< random number generator to draw from the solution pool */
	   SCIP_SOL*             selsol;             /**< best selected solution by crossover as reference point */
	};
	
	/** dins neighborhood data structure */
	struct data_dins
	{
	   int                   npoolsols;          /**< number of pool solutions where binary solution values must agree */
	};
	
	struct data_trustregion
	{
	   SCIP_Real             violpenalty;        /**< the penalty for violating the trust region */
	};
	
	/** primal heuristic data */
	struct SCIP_HeurData
	{
	   NH**                  neighborhoods;      /**< array of neighborhoods */
	   SCIP_BANDIT*          bandit;             /**< bandit algorithm */
	   SCIP_SOL*             lastcallsol;        /**< incumbent when the heuristic was last called */
	   char*                 rewardfilename;     /**< file name to store all rewards and the selection of the bandit */
	   FILE*                 rewardfile;         /**< reward file pointer, or NULL */
	   SCIP_Longint          nodesoffset;        /**< offset added to the nodes budget */
	   SCIP_Longint          maxnodes;           /**< maximum number of nodes in a single sub-SCIP */
	   SCIP_Longint          targetnodes;        /**< targeted number of nodes to start a sub-SCIP */
	   SCIP_Longint          minnodes;           /**< minimum number of nodes required to start a sub-SCIP */
	   SCIP_Longint          usednodes;          /**< total number of nodes already spent in sub-SCIPs */
	   SCIP_Longint          waitingnodes;       /**< number of nodes since last incumbent solution that the heuristic should wait */
	   SCIP_Real             nodesquot;          /**< fraction of nodes compared to the main SCIP for budget computation */
	   SCIP_Real             nodesquotmin;       /**< lower bound on fraction of nodes compared to the main SCIP for budget computation */
	   SCIP_Real             startminimprove;    /**< initial factor by which ALNS should at least improve the incumbent */
	   SCIP_Real             minimprovelow;      /**< lower threshold for the minimal improvement over the incumbent */
	   SCIP_Real             minimprovehigh;     /**< upper bound for the minimal improvement over the incumbent */
	   SCIP_Real             minimprove;         /**< factor by which ALNS should at least improve the incumbent */
	   SCIP_Real             lplimfac;           /**< limit fraction of LPs per node to interrupt sub-SCIP */
	   SCIP_Real             exp3_gamma;         /**< weight between uniform (gamma ~ 1) and weight driven (gamma ~ 0) probability distribution for exp3 */
	   SCIP_Real             exp3_beta;          /**< reward offset between 0 and 1 at every observation for exp3 */
	   SCIP_Real             epsgreedy_eps;      /**< increase exploration in epsilon-greedy bandit algorithm */
	   SCIP_Real             ucb_alpha;          /**< parameter to increase the confidence width in UCB */
	   SCIP_Real             rewardcontrol;      /**< reward control to increase the weight of the simple solution indicator
	                                               *  and decrease the weight of the closed gap reward */
	   SCIP_Real             targetnodefactor;   /**< factor by which target node number is eventually increased */
	   SCIP_Real             rewardbaseline;     /**< the reward baseline to separate successful and failed calls */
	   SCIP_Real             fixtol;             /**< tolerance by which the fixing rate may be missed without generic fixing */
	   SCIP_Real             unfixtol;           /**< tolerance by which the fixing rate may be exceeded without generic unfixing */
	   int                   nneighborhoods;     /**< number of neighborhoods */
	   int                   nactiveneighborhoods;/**< number of active neighborhoods */
	   int                   ninitneighborhoods; /**< neighborhoods that were used at least one time */
	   int                   nsolslim;           /**< limit on the number of improving solutions in a sub-SCIP call */
	   int                   seed;               /**< initial random seed for bandit algorithms and random decisions by neighborhoods */
	   int                   currneighborhood;   /**< index of currently selected neighborhood */
	   int                   ndelayedcalls;      /**< the number of delayed calls */
	   int                   maxcallssamesol;    /**< number of allowed executions of the heuristic on the same incumbent solution
	                                              *   (-1: no limit, 0: number of active neighborhoods) */
	   SCIP_Longint          firstcallthissol;   /**< counter for the number of calls on this incumbent */
	   char                  banditalgo;         /**< the bandit algorithm: (u)pper confidence bounds, (e)xp.3, epsilon (g)reedy */
	   SCIP_Bool             useredcost;         /**< should reduced cost scores be used for variable prioritization? */
	   SCIP_Bool             usedistances;       /**< should distances from fixed variables be used for variable prioritization */
	   SCIP_Bool             usepscost;          /**< should pseudo cost scores be used for variable prioritization? */
	   SCIP_Bool             domorefixings;      /**< should the ALNS heuristic do more fixings by itself based on variable prioritization
	                                               *  until the target fixing rate is reached? */
	   SCIP_Bool             adjustfixingrate;   /**< should the heuristic adjust the target fixing rate based on the success? */
	   SCIP_Bool             usesubscipheurs;    /**< should the heuristic activate other sub-SCIP heuristics during its search?  */
	   SCIP_Bool             adjustminimprove;   /**< should the factor by which the minimum improvement is bound be dynamically updated? */
	   SCIP_Bool             adjusttargetnodes;  /**< should the target nodes be dynamically adjusted? */
	   SCIP_Bool             resetweights;       /**< should the bandit algorithms be reset when a new problem is read? */
	   SCIP_Bool             subsciprandseeds;   /**< should random seeds of sub-SCIPs be altered to increase diversification? */
	   SCIP_Bool             scalebyeffort;      /**< should the reward be scaled by the effort? */
	   SCIP_Bool             copycuts;           /**< should cutting planes be copied to the sub-SCIP? */
	   SCIP_Bool             uselocalredcost;    /**< should local reduced costs be used for generic (un)fixing? */
	   SCIP_Bool             initduringroot;     /**< should the heuristic be executed multiple times during the root node? */
	};
	
	/** event handler data */
	struct SCIP_EventData
	{
	   SCIP_VAR**            subvars;            /**< the variables of the subproblem */
	   SCIP*                 sourcescip;         /**< original SCIP data structure */
	   SCIP_HEUR*            heur;               /**< alns heuristic structure */
	   SCIP_Longint          nodelimit;          /**< node limit of the run */
	   SCIP_Real             lplimfac;           /**< limit fraction of LPs per node to interrupt sub-SCIP */
	   NH_STATS*             runstats;           /**< run statistics for the current neighborhood */
	   SCIP_Bool             allrewardsmode;     /**< true if solutions should only be checked for reward comparisons */
	};
	
	/** represents limits for the sub-SCIP solving process */
	struct SolveLimits
	{
	   SCIP_Longint          nodelimit;          /**< maximum number of solving nodes for the sub-SCIP */
	   SCIP_Real             memorylimit;        /**< memory limit for the sub-SCIP */
	   SCIP_Real             timelimit;          /**< time limit for the sub-SCIP */
	   SCIP_Longint          stallnodes;         /**< maximum number of nodes without (primal) stalling */
	};
	
	typedef struct SolveLimits SOLVELIMITS;
	
	/** data structure that can be used for variable prioritization for additional fixings */
	struct VarPrio
	{
	   SCIP*                 scip;               /**< SCIP data structure */
	   SCIP_Real*            randscores;         /**< random scores for prioritization */
	   int*                  distances;          /**< breadth-first distances from already fixed variables */
	   SCIP_Real*            redcostscores;      /**< reduced cost scores for fixing a variable to a reference value */
	   SCIP_Real*            pscostscores;       /**< pseudocost scores for fixing a variable to a reference value */
	   unsigned int          useredcost:1;       /**< should reduced cost scores be used for variable prioritization? */
	   unsigned int          usedistances:1;     /**< should distances from fixed variables be used for variable prioritization */
	   unsigned int          usepscost:1;        /**< should pseudo cost scores be used for variable prioritization? */
	};
	
	/*
	 * Local methods
	 */
	
	/** Reset target fixing rate */
	static
	SCIP_RETCODE resetFixingRate(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH_FIXINGRATE*        fixingrate          /**< heuristic fixing rate */
	   )
	{
	   assert(scip != NULL);
	   assert(fixingrate != NULL);
	   fixingrate->increment = FIXINGRATE_STARTINC;
	
	   /* always start with the most conservative value */
	   fixingrate->targetfixingrate = fixingrate->maxfixingrate;
	
	   return SCIP_OKAY;
	}
	
	/** reset the currently active neighborhood */
	static
	void resetCurrentNeighborhood(
	   SCIP_HEURDATA*        heurdata
	   )
	{
	   assert(heurdata != NULL);
	   heurdata->currneighborhood = -1;
	   heurdata->ndelayedcalls = 0;
	}
	
	/** update increment for fixing rate */
	static
	void updateFixingRateIncrement(
	   NH_FIXINGRATE*        fx                  /**< fixing rate */
	   )
	{
	   fx->increment *= FIXINGRATE_DECAY;
	   fx->increment = MAX(fx->increment, LRATEMIN);
	}
	
	
	/** increase fixing rate
	 *
	 *  decrease also the rate by which the target fixing rate is adjusted
	 */
	static
	void increaseFixingRate(
	   NH_FIXINGRATE*        fx                  /**< fixing rate */
	   )
	{
	   fx->targetfixingrate += fx->increment;
	   fx->targetfixingrate = MIN(fx->targetfixingrate, fx->maxfixingrate);
	}
	
	/** decrease fixing rate
	 *
	 *  decrease also the rate by which the target fixing rate is adjusted
	 */
	static
	void decreaseFixingRate(
	   NH_FIXINGRATE*        fx                  /**< fixing rate */
	   )
	{
	   fx->targetfixingrate -= fx->increment;
	   fx->targetfixingrate = MAX(fx->targetfixingrate, fx->minfixingrate);
	}
	
	/** update fixing rate based on the results of the current run */
	static
	void updateFixingRate(
	   NH*                   neighborhood,       /**< neighborhood */
	   SCIP_STATUS           subscipstatus,      /**< status of the sub-SCIP run */
	   NH_STATS*             runstats            /**< run statistics for this run */
	   )
	{
	   NH_FIXINGRATE* fx;
	
	   fx = &neighborhood->fixingrate;
	
	   switch (subscipstatus)
	   {
	   case SCIP_STATUS_OPTIMAL:
	   case SCIP_STATUS_INFEASIBLE:
	   case SCIP_STATUS_INFORUNBD:
	   case SCIP_STATUS_SOLLIMIT:
	   case SCIP_STATUS_BESTSOLLIMIT:
	      /* decrease the fixing rate (make subproblem harder) */
	      decreaseFixingRate(fx);
	      break;
	   case SCIP_STATUS_STALLNODELIMIT:
	   case SCIP_STATUS_USERINTERRUPT:
	   case SCIP_STATUS_TERMINATE:
	   case SCIP_STATUS_NODELIMIT:
	      /* increase the fixing rate (make the subproblem easier) only if no solution was found */
	      if( runstats->nbestsolsfound <= 0 )
	         increaseFixingRate(fx);
	      break;
	      /* fall through cases to please lint */
	   case SCIP_STATUS_UNKNOWN:
	   case SCIP_STATUS_TOTALNODELIMIT:
	   case SCIP_STATUS_TIMELIMIT:
	   case SCIP_STATUS_MEMLIMIT:
	   case SCIP_STATUS_GAPLIMIT:
	   case SCIP_STATUS_RESTARTLIMIT:
	   case SCIP_STATUS_UNBOUNDED:
	   default:
	      break;
	   }
	
	   updateFixingRateIncrement(fx);
	}
	
	/** increase target node limit */
	static
	void increaseTargetNodeLimit(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
	   )
	{
	   heurdata->targetnodes = (SCIP_Longint)(heurdata->targetnodes * heurdata->targetnodefactor) + 1;
	
	   /* respect upper and lower parametrized bounds on targetnodes */
	   if( heurdata->targetnodes > heurdata->maxnodes )
	      heurdata->targetnodes = heurdata->maxnodes;
	}
	
	/** reset target node limit */
	static
	void resetTargetNodeLimit(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
	   )
	{
	   heurdata->targetnodes = heurdata->minnodes;
	}
	
	/** update target node limit based on the current run results */
	static
	void updateTargetNodeLimit(
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
	   NH_STATS*             runstats,           /**< statistics of the run */
	   SCIP_STATUS           subscipstatus       /**< status of the sub-SCIP run */
	   )
	{
	   switch (subscipstatus)
	   {
	   case SCIP_STATUS_STALLNODELIMIT:
	   case SCIP_STATUS_NODELIMIT:
	      /* the subproblem could be explored more */
	      if( runstats->nbestsolsfound == 0 )
	         increaseTargetNodeLimit(heurdata);
	      break;
	   case SCIP_STATUS_OPTIMAL:
	   case SCIP_STATUS_INFEASIBLE:
	   case SCIP_STATUS_INFORUNBD:
	   case SCIP_STATUS_SOLLIMIT:
	   case SCIP_STATUS_BESTSOLLIMIT:
	      break;
	   case SCIP_STATUS_USERINTERRUPT:
	   case SCIP_STATUS_TERMINATE:
	   case SCIP_STATUS_UNKNOWN:
	   case SCIP_STATUS_TOTALNODELIMIT:
	   case SCIP_STATUS_TIMELIMIT:
	   case SCIP_STATUS_MEMLIMIT:
	   case SCIP_STATUS_GAPLIMIT:
	   case SCIP_STATUS_RESTARTLIMIT:
	   case SCIP_STATUS_UNBOUNDED:
	      break;
	   default:
	      break;
	   }
	}
	
	/** reset the minimum improvement for the sub-SCIPs */
	static
	void resetMinimumImprovement(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
	   )
	{
	   assert(heurdata != NULL);
	   heurdata->minimprove = heurdata->startminimprove;
	}
	
	/** increase minimum improvement for the sub-SCIPs */
	static
	void increaseMinimumImprovement(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
	   )
	{
	   assert(heurdata != NULL);
	
	   heurdata->minimprove *= MINIMPROVEFAC;
	   heurdata->minimprove = MIN(heurdata->minimprove, heurdata->minimprovehigh);
	}
	
	/** decrease the minimum improvement for the sub-SCIPs */
	static
	void decreaseMinimumImprovement(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data */
	   )
	{
	   assert(heurdata != NULL);
	
	   heurdata->minimprove /= MINIMPROVEFAC;
	   SCIPdebugMessage("%.4f", heurdata->minimprovelow);
	   heurdata->minimprove = MAX(heurdata->minimprove, heurdata->minimprovelow);
	}
	
	/** update the minimum improvement based on the status of the sub-SCIP */
	static
	void updateMinimumImprovement(
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
	   SCIP_STATUS           subscipstatus,      /**< status of the sub-SCIP run */
	   NH_STATS*             runstats            /**< run statistics for this run */
	   )
	{
	   assert(heurdata != NULL);
	
	   /* if the sub-SCIP status was infeasible, we rather want to make the sub-SCIP easier
	    * with a smaller minimum improvement.
	    *
	    * If a solution limit was reached, we may, set it higher.
	    */
	   switch (subscipstatus)
	   {
	   case SCIP_STATUS_INFEASIBLE:
	   case SCIP_STATUS_INFORUNBD:
	      /* subproblem was infeasible, probably due to the minimum improvement -> decrease minimum improvement */
	      decreaseMinimumImprovement(heurdata);
	
	      break;
	   case SCIP_STATUS_SOLLIMIT:
	   case SCIP_STATUS_BESTSOLLIMIT:
	   case SCIP_STATUS_OPTIMAL:
	      /* subproblem could be optimally solved -> try higher minimum improvement */
	      increaseMinimumImprovement(heurdata);
	      break;
	   case SCIP_STATUS_NODELIMIT:
	   case SCIP_STATUS_STALLNODELIMIT:
	   case SCIP_STATUS_USERINTERRUPT:
	      /* subproblem was too hard, decrease minimum improvement */
	      if( runstats->nbestsolsfound <= 0 )
	         decreaseMinimumImprovement(heurdata);
	      break;
	   case SCIP_STATUS_UNKNOWN:
	   case SCIP_STATUS_TOTALNODELIMIT:
	   case SCIP_STATUS_TIMELIMIT:
	   case SCIP_STATUS_MEMLIMIT:
	   case SCIP_STATUS_GAPLIMIT:
	   case SCIP_STATUS_RESTARTLIMIT:
	   case SCIP_STATUS_UNBOUNDED:
	   case SCIP_STATUS_TERMINATE:
	   default:
	      break;
	   }
	}
	
	/** Reset neighborhood statistics */
	static
	SCIP_RETCODE neighborhoodStatsReset(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH_STATS*             stats               /**< neighborhood statistics */
	   )
	{
	   assert(scip != NULL);
	   assert(stats != NULL);
	
	   stats->nbestsolsfound = 0;
	   stats->nruns = 0;
	   stats->nrunsbestsol = 0;
	   stats->nsolsfound = 0;
	   stats->usednodes = 0L;
	   stats->nfixings = 0L;
	
	   BMSclearMemoryArray(stats->statushist, NHISTENTRIES);
	
	   SCIP_CALL( SCIPresetClock(scip, stats->setupclock) );
	   SCIP_CALL( SCIPresetClock(scip, stats->submipclock) );
	
	   return SCIP_OKAY;
	}
	
	/** create a neighborhood of the specified name and include it into the ALNS heuristic */
	static
	SCIP_RETCODE alnsIncludeNeighborhood(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS heuristic */
	   NH**                  neighborhood,       /**< pointer to store the neighborhood */
	   const char*           name,               /**< name for this neighborhood */
	   SCIP_Real             minfixingrate,      /**< default value for minfixingrate parameter of this neighborhood */
	   SCIP_Real             maxfixingrate,      /**< default value for maxfixingrate parameter of this neighborhood */
	   SCIP_Bool             active,             /**< default value for active parameter of this neighborhood */
	   SCIP_Real             priority,           /**< positive call priority to initialize bandit algorithms */
	   DECL_VARFIXINGS       ((*varfixings)),    /**< variable fixing callback for this neighborhood, or NULL */
	   DECL_CHANGESUBSCIP    ((*changesubscip)), /**< subscip changes callback for this neighborhood, or NULL */
	   DECL_NHINIT           ((*nhinit)),        /**< initialization callback for neighborhood, or NULL */
	   DECL_NHEXIT           ((*nhexit)),        /**< deinitialization callback for neighborhood, or NULL */
	   DECL_NHFREE           ((*nhfree)),        /**< deinitialization callback before SCIP is freed, or NULL */
	   DECL_NHREFSOL         ((*nhrefsol)),      /**< callback function to return a reference solution for further fixings, or NULL */
	   DECL_NHDEACTIVATE     ((*nhdeactivate))   /**< callback function to deactivate neighborhoods on problems where they are irrelevant, or NULL if neighborhood is always active */
	   )
	{
	   char paramname[SCIP_MAXSTRLEN];
	
	   assert(scip != NULL);
	   assert(heurdata != NULL);
	   assert(neighborhood != NULL);
	   assert(name != NULL);
	
	   SCIP_CALL( SCIPallocBlockMemory(scip, neighborhood) );
	   assert(*neighborhood != NULL);
	
	   SCIP_ALLOC( BMSduplicateMemoryArray(&(*neighborhood)->name, name, strlen(name)+1) );
	
	   SCIP_CALL( SCIPcreateClock(scip, &(*neighborhood)->stats.setupclock) );
	   SCIP_CALL( SCIPcreateClock(scip, &(*neighborhood)->stats.submipclock) );
	
	   (*neighborhood)->changesubscip = changesubscip;
	   (*neighborhood)->varfixings = varfixings;
	   (*neighborhood)->nhinit = nhinit;
	   (*neighborhood)->nhexit = nhexit;
	   (*neighborhood)->nhfree = nhfree;
	   (*neighborhood)->nhrefsol = nhrefsol;
	   (*neighborhood)->nhdeactivate = nhdeactivate;
	
	   /* add parameters for this neighborhood */
	   (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "heuristics/alns/%s/minfixingrate", name);
	   SCIP_CALL( SCIPaddRealParam(scip, paramname, "minimum fixing rate for this neighborhood",
	         &(*neighborhood)->fixingrate.minfixingrate, TRUE, minfixingrate, 0.0, 1.0, NULL, NULL) );
	   (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "heuristics/alns/%s/maxfixingrate", name);
	   SCIP_CALL( SCIPaddRealParam(scip, paramname, "maximum fixing rate for this neighborhood",
	         &(*neighborhood)->fixingrate.maxfixingrate, TRUE, maxfixingrate, 0.0, 1.0, NULL, NULL) );
	   (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "heuristics/alns/%s/active", name);
	   SCIP_CALL( SCIPaddBoolParam(scip, paramname, "is this neighborhood active?",
	         &(*neighborhood)->active, TRUE, active, NULL, NULL) );
	   (void) SCIPsnprintf(paramname, SCIP_MAXSTRLEN, "heuristics/alns/%s/priority", name);
	   SCIP_CALL( SCIPaddRealParam(scip, paramname, "positive call priority to initialize bandit algorithms",
	         &(*neighborhood)->priority, TRUE, priority, 1e-2, 1.0, NULL, NULL) );
	
	   /* add the neighborhood to the ALNS heuristic */
	   heurdata->neighborhoods[heurdata->nneighborhoods++] = (*neighborhood);
	
	   return SCIP_OKAY;
	}
	
	/** release all data and free neighborhood */
	static
	SCIP_RETCODE alnsFreeNeighborhood(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH**                  neighborhood        /**< pointer to neighborhood that should be freed */
	   )
	{
	   NH* nhptr;
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	
	   nhptr = *neighborhood;
	   assert(nhptr != NULL);
	
	   BMSfreeMemoryArray(&nhptr->name);
	
	   /* release further, neighborhood specific data structures */
	   if( nhptr->nhfree != NULL )
	   {
	      SCIP_CALL( nhptr->nhfree(scip, nhptr) );
	   }
	
	   SCIP_CALL( SCIPfreeClock(scip, &nhptr->stats.setupclock) );
	   SCIP_CALL( SCIPfreeClock(scip, &nhptr->stats.submipclock) );
	
	   SCIPfreeBlockMemory(scip, neighborhood);
	   *neighborhood = NULL;
	
	   return SCIP_OKAY;
	}
	
	/** initialize neighborhood specific data */
	static
	SCIP_RETCODE neighborhoodInit(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH*                   neighborhood        /**< neighborhood to initialize */
	   )
	{
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	
	   /* call the init callback of the neighborhood */
	   if( neighborhood->nhinit != NULL )
	   {
	      SCIP_CALL( neighborhood->nhinit(scip, neighborhood) );
	   }
	
	   return SCIP_OKAY;
	}
	
	/** deinitialize neighborhood specific data */
	static
	SCIP_RETCODE neighborhoodExit(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH*                   neighborhood        /**< neighborhood to initialize */
	   )
	{
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	
	   if( neighborhood->nhexit != NULL )
	   {
	      SCIP_CALL( neighborhood->nhexit(scip, neighborhood) );
	   }
	
	   return SCIP_OKAY;
	}
	
	/** creates a new solution for the original problem by copying the solution of the subproblem */
	static
	SCIP_RETCODE transferSolution(
	   SCIP*                 subscip,            /**< SCIP data structure of the subproblem */
	   SCIP_EVENTDATA*       eventdata           /**< event handler data */
	   )
	{
	   SCIP*      sourcescip;         /* original SCIP data structure */
	   SCIP_VAR** subvars;            /* the variables of the subproblem */
	   SCIP_HEUR* heur;               /* alns heuristic structure */
	   SCIP_SOL*  subsol;             /* solution of the subproblem */
	   SCIP_SOL*  newsol;             /* solution to be created for the original problem */
	   SCIP_Bool  success;
	   NH_STATS*  runstats;
	   SCIP_SOL*  oldbestsol;
	
	   assert(subscip != NULL);
	
	   subsol = SCIPgetBestSol(subscip);
	   assert(subsol != NULL);
	
	   sourcescip = eventdata->sourcescip;
	   subvars = eventdata->subvars;
	   heur = eventdata->heur;
	   runstats = eventdata->runstats;
	   assert(sourcescip != NULL);
	   assert(sourcescip != subscip);
	   assert(heur != NULL);
	   assert(subvars != NULL);
	   assert(runstats != NULL);
	
	   SCIP_CALL( SCIPtranslateSubSol(sourcescip, subscip, subsol, heur, subvars, &newsol) );
	
	   oldbestsol = SCIPgetBestSol(sourcescip);
	
	   /* in the special, experimental all rewards mode, the solution is only checked for feasibility
	    * but not stored
	    */
	   if( eventdata->allrewardsmode )
	   {
	      SCIP_CALL( SCIPcheckSol(sourcescip, newsol, FALSE, FALSE, TRUE, TRUE, TRUE, &success) );
	
	      if( success )
	      {
	         runstats->nsolsfound++;
	         if( SCIPgetSolTransObj(sourcescip, newsol) < SCIPgetCutoffbound(sourcescip) )
	            runstats->nbestsolsfound++;
	      }
	
	      SCIP_CALL( SCIPfreeSol(sourcescip, &newsol) );
	   }
	   else
	   {
	      /* try to add new solution to scip and free it immediately */
	      SCIP_CALL( SCIPtrySolFree(sourcescip, &newsol, FALSE, FALSE, TRUE, TRUE, TRUE, &success) );
	
	      if( success )
	      {
	         runstats->nsolsfound++;
	         if( SCIPgetBestSol(sourcescip) != oldbestsol )
	            runstats->nbestsolsfound++;
	      }
	   }
	
	   /* update new upper bound for reward later */
	   runstats->newupperbound = SCIPgetUpperbound(sourcescip);
	
	   return SCIP_OKAY;
	}
	
	
	/* ---------------- Callback methods of event handler ---------------- */
	
	/** execution callback of the event handler
	 *
	 * transfer new solutions or interrupt the solving process manually
	 */
	static
	SCIP_DECL_EVENTEXEC(eventExecAlns)
	{
	   assert(eventhdlr != NULL);
	   assert(eventdata != NULL);
	   assert(strcmp(SCIPeventhdlrGetName(eventhdlr), EVENTHDLR_NAME) == 0);
	   assert(event != NULL);
	   assert(SCIPeventGetType(event) & SCIP_EVENTTYPE_ALNS);
	   assert(eventdata != NULL);
	
	   /* treat the different atomic events */
	   switch( SCIPeventGetType(event) )
	   {
	   case SCIP_EVENTTYPE_SOLFOUND:
	   case SCIP_EVENTTYPE_BESTSOLFOUND:
	      /* try to transfer the solution to the original SCIP */
	      SCIP_CALL( transferSolution(scip, eventdata) );
	      break;
	   case SCIP_EVENTTYPE_LPSOLVED:
	      /* interrupt solution process of sub-SCIP */
	      if( SCIPgetNLPs(scip) > eventdata->lplimfac * eventdata->nodelimit )
	      {
	         SCIPdebugMsg(scip, "interrupt after  %" SCIP_LONGINT_FORMAT " LPs\n", SCIPgetNLPs(scip));
	         SCIP_CALL( SCIPinterruptSolve(scip) );
	      }
	      break;
	   default:
	      break;
	   }
	
	   return SCIP_OKAY;
	}
	
	/** initialize neighborhood statistics before the next run */
	static
	void initRunStats(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH_STATS*             stats               /**< run statistics */
	   )
	{
	   stats->nbestsolsfound = 0;
	   stats->nsolsfound = 0;
	   stats->usednodes = 0L;
	   stats->nfixings = 0;
	   stats->oldupperbound = SCIPgetUpperbound(scip);
	   stats->newupperbound = SCIPgetUpperbound(scip);
	}
	
	/** update run stats after the sub SCIP was solved */
	static
	void updateRunStats(
	   NH_STATS*             stats,              /**< run statistics */
	   SCIP*                 subscip             /**< sub-SCIP instance, or NULL */
	   )
	{
	   /* treat an untransformed subscip as if none was created */
	   if( subscip != NULL && ! SCIPisTransformed(subscip) )
	      subscip = NULL;
	
	   stats->usednodes = subscip != NULL ? SCIPgetNNodes(subscip) : 0L;
	}
	
	/** get the histogram index for this status */
	static
	int getHistIndex(
	   SCIP_STATUS           subscipstatus       /**< sub-SCIP status */
	   )
	{
	   switch (subscipstatus)
	   {
	   case SCIP_STATUS_OPTIMAL:
	      return (int)HIDX_OPT;
	   case SCIP_STATUS_INFEASIBLE:
	      return (int)HIDX_INFEAS;
	   case SCIP_STATUS_NODELIMIT:
	      return (int)HIDX_NODELIM;
	   case SCIP_STATUS_STALLNODELIMIT:
	      return (int)HIDX_STALLNODE;
	   case SCIP_STATUS_SOLLIMIT:
	   case SCIP_STATUS_BESTSOLLIMIT:
	      return (int)HIDX_SOLLIM;
	   case SCIP_STATUS_USERINTERRUPT:
	      return (int)HIDX_USR;
	   default:
	      return (int)HIDX_OTHER;
	   } /*lint !e788*/
	}
	
	/** print neighborhood statistics */
	static
	void printNeighborhoodStatistics(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data */
	   FILE*                 file                /**< file handle, or NULL for standard out */
	   )
	{
	   int i;
	   int j;
	   HISTINDEX statusses[] = {HIDX_OPT, HIDX_INFEAS, HIDX_NODELIM, HIDX_STALLNODE, HIDX_SOLLIM, HIDX_USR, HIDX_OTHER};
	
	   SCIPinfoMessage(scip, file, "Neighborhoods      : %10s %10s %10s %10s %10s %10s %10s %10s %10s %10s %4s %4s %4s %4s %4s %4s %4s %4s\n",
	            "Calls", "SetupTime", "SolveTime", "SolveNodes", "Sols", "Best", "Exp3", "EpsGreedy", "UCB", "TgtFixRate",
	            "Opt", "Inf", "Node", "Stal", "Sol", "Usr", "Othr", "Actv");
	
	   /* loop over neighborhoods and fill in statistics */
	   for( i = 0; i < heurdata->nneighborhoods; ++i )
	   {
	      NH* neighborhood;
	      SCIP_Real proba;
	      SCIP_Real ucb;
	      SCIP_Real epsgreedyweight;
	
	      neighborhood = heurdata->neighborhoods[i];
	      SCIPinfoMessage(scip, file, "  %-17s:", neighborhood->name);
	      SCIPinfoMessage(scip, file, " %10d", neighborhood->stats.nruns);
	      SCIPinfoMessage(scip, file, " %10.2f", SCIPgetClockTime(scip, neighborhood->stats.setupclock) );
	      SCIPinfoMessage(scip, file, " %10.2f", SCIPgetClockTime(scip, neighborhood->stats.submipclock) );
	      SCIPinfoMessage(scip, file, " %10" SCIP_LONGINT_FORMAT, neighborhood->stats.usednodes );
	      SCIPinfoMessage(scip, file, " %10" SCIP_LONGINT_FORMAT, neighborhood->stats.nsolsfound);
	      SCIPinfoMessage(scip, file, " %10" SCIP_LONGINT_FORMAT, neighborhood->stats.nbestsolsfound);
	
	      proba = 0.0;
	      ucb = 1.0;
	      epsgreedyweight = -1.0;
	
	      if( heurdata->bandit != NULL && i < heurdata->nactiveneighborhoods )
	      {
	         switch (heurdata->banditalgo)
	         {
	         case 'u':
	            ucb = SCIPgetConfidenceBoundUcb(heurdata->bandit, i);
	            break;
	         case 'g':
	            epsgreedyweight = SCIPgetWeightsEpsgreedy(heurdata->bandit)[i];
	            break;
	         case 'e':
	            proba = SCIPgetProbabilityExp3(heurdata->bandit, i);
	            break;
	         default:
	            break;
	         }
	      }
	
	      SCIPinfoMessage(scip, file, " %10.5f", proba);
	      SCIPinfoMessage(scip, file, " %10.5f", epsgreedyweight);
	      SCIPinfoMessage(scip, file, " %10.5f", ucb);
	      SCIPinfoMessage(scip, file, " %10.3f", neighborhood->fixingrate.targetfixingrate);
	
	      /* loop over status histogram */
	      for( j = 0; j < NHISTENTRIES; ++j )
	         SCIPinfoMessage(scip, file, " %4d", neighborhood->stats.statushist[statusses[j]]);
	
	      SCIPinfoMessage(scip, file, " %4d", i < heurdata->nactiveneighborhoods ? 1 : 0);
	      SCIPinfoMessage(scip, file, "\n");
	   }
	}
	
	/** update the statistics of the neighborhood based on the sub-SCIP run */
	static
	void updateNeighborhoodStats(
	   NH_STATS*             runstats,           /**< run statistics */
	   NH*                   neighborhood,       /**< the selected neighborhood */
	   SCIP_STATUS           subscipstatus       /**< status of the sub-SCIP solve */
	   )
	{  /*lint --e{715}*/
	   NH_STATS* stats;
	   stats = &neighborhood->stats;
	
	   /* copy run statistics into neighborhood statistics */
	   stats->nbestsolsfound += runstats->nbestsolsfound;
	   stats->nsolsfound += runstats->nsolsfound;
	   stats->usednodes += runstats->usednodes;
	   stats->nruns += 1;
	
	   if( runstats->nbestsolsfound > 0 )
	      stats->nrunsbestsol += DEFAULT_BESTSOLWEIGHT;
	   else if( runstats->nsolsfound > 0 )
	      stats->nrunsbestsol++;
	
	   /* update the counter for the subscip status */
	   ++stats->statushist[getHistIndex(subscipstatus)];
	}
	
	/** sort callback for variable pointers using the ALNS variable prioritization
	 *
	 *  the variable prioritization works hierarchically as follows. A variable
	 *  a has the higher priority over b iff
	 *
	 *  - variable distances should be used and a has a smaller distance than b
	 *  - variable reduced costs should be used and a has a smaller score than b
	 *  - variable pseudo costs should be used and a has a smaller score than b
	 *  - based on previously assigned random scores
	 *
	 *  @note: distances are context-based. For fixing more variables,
	 *  distances are initialized from the already fixed variables.
	 *  For unfixing variables, distances are initialized starting
	 *  from the unfixed variables
	 */
	static
	SCIP_DECL_SORTINDCOMP(sortIndCompAlns)
	{  /*lint --e{715}*/
	   VARPRIO* varprio;
	
	   varprio = (VARPRIO*)dataptr;
	   assert(varprio != NULL);
	   assert(varprio->randscores != NULL);
	
	   if( ind1 == ind2 )
	      return 0;
	
	   /* priority is on distances, if enabled. The variable which is closer in a breadth-first search sense to
	    * the already fixed variables has precedence */
	   if( varprio->usedistances )
	   {
	      int dist1;
	      int dist2;
	
	      dist1 = varprio->distances[ind1];
	      dist2 = varprio->distances[ind2];
	
	      if( dist1 < 0 )
	         dist1 = INT_MAX;
	
	      if( dist2 < 0 )
	         dist2 = INT_MAX;
	
	      assert(varprio->distances != NULL);
	      if( dist1 < dist2 )
	         return -1;
	      else if( dist1 > dist2 )
	         return 1;
	   }
	
	   assert(! varprio->usedistances || varprio->distances[ind1] == varprio->distances[ind2]);
	
	   /* if the indices tie considering distances or distances are disabled -> use reduced cost information instead */
	   if( varprio->useredcost )
	   {
	      assert(varprio->redcostscores != NULL);
	
	      if( varprio->redcostscores[ind1] < varprio->redcostscores[ind2] )
	         return -1;
	      else if( varprio->redcostscores[ind1] > varprio->redcostscores[ind2] )
	         return 1;
	   }
	
	   /* use pseudo cost scores if reduced costs are disabled or a tie was found */
	   if( varprio->usepscost )
	   {
	      assert(varprio->pscostscores != NULL);
	
	      /* prefer the variable with smaller pseudocost score */
	      if( varprio->pscostscores[ind1] < varprio->pscostscores[ind2] )
	         return -1;
	      else if( varprio->pscostscores[ind1] > varprio->pscostscores[ind2] )
	         return 1;
	   }
	
	   if( varprio->randscores[ind1] < varprio->randscores[ind2] )
	      return -1;
	   else if( varprio->randscores[ind1] > varprio->randscores[ind2] )
	      return 1;
	
	   return ind1 - ind2;
	}
	
	/** Compute the reduced cost score for this variable in the reference solution */
	static
	SCIP_Real getVariableRedcostScore(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_VAR*             var,                /**< the variable for which the score should be computed */
	   SCIP_Real             refsolval,          /**< solution value in reference solution */
	   SCIP_Bool             uselocalredcost     /**< should local reduced costs be used for generic (un)fixing? */
	   )
	{
	   SCIP_Real bestbound;
	   SCIP_Real redcost;
	   SCIP_Real score;
	   assert(scip != NULL);
	   assert(var != NULL);
	
	   /* prefer column variables */
	   if( SCIPvarGetStatus(var) != SCIP_VARSTATUS_COLUMN )
	      return SCIPinfinity(scip);
	
	   if( ! uselocalredcost )
	   {
	      redcost = SCIPvarGetBestRootRedcost(var);
	
	      bestbound = SCIPvarGetBestRootSol(var);
	
	      /* using global reduced costs, the two factors yield a nonnegative score within tolerances */
	      assert(SCIPisDualfeasZero(scip, redcost)
	         || (SCIPisDualfeasNegative(scip, redcost) && ! SCIPisFeasPositive(scip, refsolval - bestbound))
	         || (SCIPisDualfeasPositive(scip, redcost) && ! SCIPisFeasNegative(scip, refsolval - bestbound)));
	   }
	   else
	   {
	      /* this can be safely asserted here, since the heuristic would not reach this point, otherwise */
	      assert(SCIPhasCurrentNodeLP(scip));
	      assert(SCIPgetLPSolstat(scip) == SCIP_LPSOLSTAT_OPTIMAL);
	
	      redcost = SCIPgetVarRedcost(scip, var);
	
	      bestbound = SCIPvarGetLPSol(var);
	   }
	
	   assert(! SCIPisInfinity(scip, REALABS(bestbound)));
	   assert(SCIPisDualfeasZero(scip, redcost) || SCIPisFeasIntegral(scip, bestbound));
	
	   score = redcost * (refsolval - bestbound);
	
	   /* max out numerical inaccuracies from global scores */
	   if( ! uselocalredcost )
	      score = MAX(score, 0.0);
	
	   return score;
	}
	
	/** get the pseudo cost score of this variable with respect to the reference solution */
	static
	SCIP_Real getVariablePscostScore(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_VAR*             var,                /**< the variable for which the score should be computed */
	   SCIP_Real             refsolval,          /**< solution value in reference solution */
	   SCIP_Bool             uselocallpsol       /**< should local LP solution be used? */
	   )
	{
	   SCIP_Real lpsolval;
	
	   assert(scip != NULL);
	   assert(var != NULL);
	
	   /* variables that aren't LP columns have no pseudocost score */
	   if( SCIPvarGetStatus(var) != SCIP_VARSTATUS_COLUMN )
	      return 0.0;
	
	   lpsolval = uselocallpsol ? SCIPvarGetLPSol(var) : SCIPvarGetRootSol(var);
	
	   /* the score is 0.0 if the values are equal */
	   if( SCIPisEQ(scip, lpsolval, refsolval) )
	      return 0.0;
	   else
	      return SCIPgetVarPseudocostVal(scip, var, refsolval - lpsolval);
	}
	
	/** add variable and solution value to buffer data structure for variable fixings. The method checks if
	 *  the value still lies within the variable bounds. The value stays unfixed otherwise.
	 */
	static
	void tryAdd2variableBuffer(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_VAR*             var,                /**< (source) SCIP variable that should be added to the buffer */
	   SCIP_Real             val,                /**< fixing value for this variable */
	   SCIP_VAR**            varbuf,             /**< variable buffer to store variables that should be fixed */
	   SCIP_Real*            valbuf,             /**< value buffer to store fixing values */
	   int*                  nfixings,           /**< pointer to number of fixed buffer variables, will be increased by 1 */
	   SCIP_Bool             integer             /**< is this an integer variable? */
	   )
	{
	   assert(SCIPisFeasIntegral(scip, val) || ! SCIPvarIsIntegral(var));
	   assert(*nfixings < SCIPgetNVars(scip));
	
	   /* round the value to its nearest integer */
	   if( integer )
	      val = SCIPfloor(scip, val + 0.5);
	
	   /* only add fixing if it is still valid within the global variable bounds. Invalidity
	    * of this solution value may come from a dual reduction that was performed after the solution from which
	    * this value originated was found
	    */
	   if( SCIPvarGetLbGlobal(var) <= val && val <= SCIPvarGetUbGlobal(var) )
	   {
	      varbuf[*nfixings] = var;
	      valbuf[*nfixings] = val;
	      ++(*nfixings);
	   }
	}
	
	/** query neighborhood for a reference solution for further fixings */
	static
	SCIP_RETCODE neighborhoodGetRefsol(
	   SCIP*                 scip,               /**< SCIP data structure */
	   NH*                   neighborhood,       /**< ALNS neighborhood data structure */
	   SCIP_SOL**            solptr              /**< solution pointer */
	   )
	{
	   assert(solptr != NULL);
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	
	   *solptr = NULL;
	   if( neighborhood->nhrefsol != NULL )
	   {
	      SCIP_RESULT result;
	      SCIP_CALL( neighborhood->nhrefsol(scip, neighborhood, solptr, &result) );
	
	      if( result == SCIP_DIDNOTFIND )
	         *solptr = NULL;
	      else
	         assert(*solptr != NULL);
	   }
	
	   return SCIP_OKAY;
	}
	
	/** fix additional variables found in feasible reference solution if the ones that the neighborhood found were not enough
	 *
	 *  use not always the best solution for the values, but a reference solution provided by the neighborhood itself
	 *
	 *  @note it may happen that the target fixing rate is not completely reached. This is the case if intermediate,
	 *  dual reductions render the solution values of the reference solution infeasible for
	 *  the current, global variable bounds.
	 */
	static
	SCIP_RETCODE alnsFixMoreVariables(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   SCIP_SOL*             refsol,             /**< feasible reference solution for more variable fixings */
	   SCIP_VAR**            varbuf,             /**< buffer array to store variables to fix */
	   SCIP_Real*            valbuf,             /**< buffer array to store fixing values */
	   int*                  nfixings,           /**< pointer to store the number of fixings */
	   int                   ntargetfixings,     /**< number of required target fixings */
	   SCIP_Bool*            success             /**< pointer to store whether the target fixings have been successfully reached */
	   )
	{
	   VARPRIO varprio;
	   SCIP_VAR** vars;
	   SCIP_Real* redcostscores;
	   SCIP_Real* pscostscores;
	   SCIP_Real* solvals;
	   SCIP_RANDNUMGEN* rng;
	   SCIP_VAR** unfixedvars;
	   SCIP_Bool* isfixed;
	   int* distances;
	   int* perm;
	   SCIP_Real* randscores;
	   int nbinvars;
	   int nintvars;
	   int nbinintvars;
	   int nvars;
	   int b;
	   int nvarstoadd;
	   int nunfixedvars;
	
	   assert(scip != NULL);
	   assert(varbuf != NULL);
	   assert(nfixings != NULL);
	   assert(success != NULL);
	   assert(heurdata != NULL);
	   assert(refsol != NULL);
	
	   *success = FALSE;
	
	   /* if the user parameter forbids more fixings, return immediately */
	   if( ! heurdata->domorefixings )
	      return SCIP_OKAY;
	
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, &nbinvars, &nintvars, NULL, NULL) );
	
	   nbinintvars = nbinvars + nintvars;
	
	   if( ntargetfixings >= nbinintvars )
	      return SCIP_OKAY;
	
	   /* determine the number of required additional fixings */
	   nvarstoadd = ntargetfixings - *nfixings;
	   if( nvarstoadd == 0 )
	      return SCIP_OKAY;
	
	   varprio.usedistances = heurdata->usedistances && (*nfixings >= 1);
	   varprio.useredcost = heurdata->useredcost;
	   varprio.usepscost = heurdata->usepscost;
	   varprio.scip = scip;
	   rng = SCIPbanditGetRandnumgen(heurdata->bandit);
	   assert(rng != NULL);
	
	   SCIP_CALL( SCIPallocBufferArray(scip, &randscores, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &perm, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &distances, nvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &redcostscores, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &solvals, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &isfixed, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &unfixedvars, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &pscostscores, nbinintvars) );
	
	   /* initialize variable graph distances from already fixed variables */
	   if( varprio.usedistances )
	   {
	      SCIP_CALL( SCIPvariablegraphBreadthFirst(scip, NULL, varbuf, *nfixings, distances, INT_MAX, INT_MAX, ntargetfixings) );
	   }
	   else
	   {
	      /* initialize all equal distances to make them irrelevant */
	      BMSclearMemoryArray(distances, nbinintvars);
	   }
	
	   BMSclearMemoryArray(isfixed, nbinintvars);
	
	   /* mark binary and integer variables if they are fixed */
	   for( b = 0; b < *nfixings; ++b )
	   {
	      int probindex;
	
	      assert(varbuf[b] != NULL);
	      probindex = SCIPvarGetProbindex(varbuf[b]);
	      assert(probindex >= 0);
	
	      if( probindex < nbinintvars )
	         isfixed[probindex] = TRUE;
	   }
	
	   SCIP_CALL( SCIPgetSolVals(scip, refsol, nbinintvars, vars, solvals) );
	
	   /* assign scores to unfixed every discrete variable of the problem */
	   nunfixedvars = 0;
	   for( b = 0; b < nbinintvars; ++b )
	   {
	      SCIP_VAR* var = vars[b];
	
	      /* filter fixed variables */
	      if( isfixed[b] )
	         continue;
	
	      /* filter variables with a solution value outside its global bounds */
	      if( solvals[b] < SCIPvarGetLbGlobal(var) - 0.5 || solvals[b] > SCIPvarGetUbGlobal(var) + 0.5 )
	         continue;
	
	      redcostscores[nunfixedvars] = getVariableRedcostScore(scip, var, solvals[b], heurdata->uselocalredcost);
	      pscostscores[nunfixedvars] = getVariablePscostScore(scip, var, solvals[b], heurdata->uselocalredcost);
	
	      unfixedvars[nunfixedvars] = var;
	      perm[nunfixedvars] = nunfixedvars;
	      randscores[nunfixedvars] = SCIPrandomGetReal(rng, 0.0, 1.0);
	
	      /* these assignments are based on the fact that nunfixedvars <= b */
	      solvals[nunfixedvars] = solvals[b];
	      distances[nunfixedvars] = distances[b];
	
	      SCIPdebugMsg(scip, "Var <%s> scores: dist %3d, red cost %15.9g, pscost %15.9g rand %6.4f\n",
	         SCIPvarGetName(var), distances[nunfixedvars], redcostscores[nunfixedvars],
	         pscostscores[nunfixedvars], randscores[nunfixedvars]);
	
	      nunfixedvars++;
	   }
	
	   /* use selection algorithm (order of the variables does not matter) for quickly completing the fixing */
	   varprio.randscores = randscores;
	   varprio.distances = distances;
	   varprio.redcostscores = redcostscores;
	   varprio.pscostscores = pscostscores;
	
	   nvarstoadd = MIN(nunfixedvars, nvarstoadd);
	
	   /* select the first nvarstoadd many variables according to the score */
	   if( nvarstoadd < nunfixedvars )
	      SCIPselectInd(perm, sortIndCompAlns, &varprio, nvarstoadd, nunfixedvars);
	
	   /* loop over the first elements of the selection defined in permutation. They represent the best variables */
	   for( b = 0; b < nvarstoadd; ++b )
	   {
	      int permindex = perm[b];
	      assert(permindex >= 0);
	      assert(permindex < nunfixedvars);
	
	      tryAdd2variableBuffer(scip, unfixedvars[permindex], solvals[permindex], varbuf, valbuf, nfixings, TRUE);
	   }
	
	   *success = TRUE;
	
	   /* free buffer arrays */
	   SCIPfreeBufferArray(scip, &pscostscores);
	   SCIPfreeBufferArray(scip, &unfixedvars);
	   SCIPfreeBufferArray(scip, &isfixed);
	   SCIPfreeBufferArray(scip, &solvals);
	   SCIPfreeBufferArray(scip, &redcostscores);
	   SCIPfreeBufferArray(scip, &distances);
	   SCIPfreeBufferArray(scip, &perm);
	   SCIPfreeBufferArray(scip, &randscores);
	
	   return SCIP_OKAY;
	}
	
	/** create the bandit algorithm for the heuristic depending on the user parameter */
	static
	SCIP_RETCODE createBandit(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data structure */
	   SCIP_Real*            priorities,         /**< call priorities for active neighborhoods */
	   unsigned int          initseed            /**< initial random seed */
	   )
	{
	   switch (heurdata->banditalgo)
	   {
	   case 'u':
	      SCIP_CALL( SCIPcreateBanditUcb(scip, &heurdata->bandit, priorities,
	            heurdata->ucb_alpha, heurdata->nactiveneighborhoods, initseed) );
	      break;
	
	   case 'e':
	      SCIP_CALL( SCIPcreateBanditExp3(scip, &heurdata->bandit, priorities,
	            heurdata->exp3_gamma, heurdata->exp3_beta, heurdata->nactiveneighborhoods, initseed) );
	      break;
	
	   case 'g':
	      SCIP_CALL( SCIPcreateBanditEpsgreedy(scip, &heurdata->bandit, priorities,
	            heurdata->epsgreedy_eps, FALSE, 0.9, 0, heurdata->nactiveneighborhoods, initseed) );
	      break;
	
	   default:
	      SCIPerrorMessage("Unknown bandit parameter %c\n", heurdata->banditalgo);
	      return SCIP_INVALIDDATA;
	   }
	
	   return SCIP_OKAY;
	}
	
	/*
	 * Callback methods of primal heuristic
	 */
	
	/** copy method for primal heuristic plugins (called when SCIP copies plugins) */
	static
	SCIP_DECL_HEURCOPY(heurCopyAlns)
	{  /*lint --e{715}*/
	   assert(scip != NULL);
	   assert(heur != NULL);
	   assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
	
	   /* call inclusion method of primal heuristic */
	   SCIP_CALL( SCIPincludeHeurAlns(scip) );
	
	   return SCIP_OKAY;
	}
	
	/** unfix some of the variables because there are too many fixed
	 *
	 *  a variable is ideally unfixed if it is close to other unfixed variables
	 *  and fixing it has a high reduced cost impact
	 */
	static
	SCIP_RETCODE alnsUnfixVariables(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   SCIP_VAR**            varbuf,             /**< buffer array to store variables to fix */
	   SCIP_Real*            valbuf,             /**< buffer array to store fixing values */
	   int*                  nfixings,           /**< pointer to store the number of fixings */
	   int                   ntargetfixings,     /**< number of required target fixings */
	   SCIP_Bool*            success             /**< pointer to store whether the target fixings have been successfully reached */
	   )
	{
	   VARPRIO varprio;
	   SCIP_Real* redcostscores;
	   SCIP_Real* pscostscores;
	   SCIP_Real* randscores;
	   SCIP_VAR** unfixedvars;
	   SCIP_VAR** varbufcpy;
	   SCIP_Real* valbufcpy;
	   SCIP_Bool* isfixedvar;
	   SCIP_VAR** vars;
	   SCIP_RANDNUMGEN* rng;
	   int* distances;
	   int* fixeddistances;
	   int* perm;
	   int nvars;
	   int i;
	   int nbinintvars;
	   int nunfixed;
	
	   *success = FALSE;
	
	   nbinintvars = SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip);
	   if( nbinintvars == 0 )
	      return SCIP_OKAY;
	
	   assert(*nfixings > 0);
	
	   nvars = SCIPgetNVars(scip);
	   SCIP_CALL( SCIPallocBufferArray(scip, &isfixedvar, nvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &unfixedvars, nbinintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &distances, nvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &fixeddistances, *nfixings) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &redcostscores, *nfixings) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &randscores, *nfixings) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &perm, *nfixings) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &pscostscores, *nfixings) );
	
	   SCIP_CALL( SCIPduplicateBufferArray(scip, &varbufcpy, varbuf, *nfixings) );
	   SCIP_CALL( SCIPduplicateBufferArray(scip, &valbufcpy, valbuf, *nfixings) );
	
	   /*
	    * collect the unfixed binary and integer variables
	    */
	   BMSclearMemoryArray(isfixedvar, nvars);
	   /* loop over fixed variables and mark their respective positions as fixed */
	   for( i = 0; i < *nfixings; ++i )
	   {
	      int probindex = SCIPvarGetProbindex(varbuf[i]);
	
	      assert(probindex >= 0);
	
	      isfixedvar[probindex] = TRUE;
	   }
	
	   nunfixed = 0;
	   vars = SCIPgetVars(scip);
	   /* collect unfixed binary and integer variables */
	   for( i = 0; i < nbinintvars; ++i )
	   {
	      if( ! isfixedvar[i] )
	         unfixedvars[nunfixed++] = vars[i];
	   }
	
	   varprio.usedistances = heurdata->usedistances && nunfixed > 0;
	
	   /* collect distances of all fixed variables from those that are not fixed */
	   if( varprio.usedistances )
	   {
	      SCIP_CALL( SCIPvariablegraphBreadthFirst(scip, NULL, unfixedvars, nunfixed, distances, INT_MAX, INT_MAX, INT_MAX) );
	
	      for( i = 0; i < *nfixings; ++i )
	      {
	         int probindex = SCIPvarGetProbindex(varbuf[i]);
	         if( probindex >= 0 )
	            fixeddistances[i] = distances[probindex];
	      }
	   }
	   else
	   {
	      BMSclearMemoryArray(fixeddistances, *nfixings);
	   }
	
	   /* collect reduced cost scores of the fixings and assign random scores */
	   rng = SCIPbanditGetRandnumgen(heurdata->bandit);
	   for( i = 0; i < *nfixings; ++i )
	   {
	      SCIP_VAR* fixedvar = varbuf[i];
	      SCIP_Real fixval = valbuf[i];
	
	      /* use negative reduced cost and pseudo cost scores to prefer variable fixings with small score */
	      redcostscores[i] = - getVariableRedcostScore(scip, fixedvar, fixval, heurdata->uselocalredcost);
	      pscostscores[i] = - getVariablePscostScore(scip, fixedvar, fixval, heurdata->uselocalredcost);
	      randscores[i] = SCIPrandomGetReal(rng, 0.0, 1.0);
	      perm[i] = i;
	
	      SCIPdebugMsg(scip, "Var <%s> scores: dist %3d, red cost %15.9g, pscost %15.9g rand %6.4f\n",
	            SCIPvarGetName(fixedvar), fixeddistances[i], redcostscores[i], pscostscores[i], randscores[i]);
	   }
	
	   varprio.distances = fixeddistances;
	   varprio.randscores = randscores;
	   varprio.redcostscores = redcostscores;
	   varprio.pscostscores = pscostscores;
	   varprio.useredcost = heurdata->useredcost;
	   varprio.usepscost = heurdata->usepscost;
	   varprio.scip = scip;
	
	   /* scores are assigned in such a way that variables with a smaller score should be fixed last */
	   SCIPselectDownInd(perm, sortIndCompAlns, &varprio, ntargetfixings, *nfixings);
	
	   /* bring the desired variables to the front of the array */
	   for( i = 0; i < ntargetfixings; ++i )
	   {
	      valbuf[i] = valbufcpy[perm[i]];
	      varbuf[i] = varbufcpy[perm[i]];
	   }
	
	   *nfixings = ntargetfixings;
	
	   /* free the buffer arrays in reverse order of allocation */
	   SCIPfreeBufferArray(scip, &valbufcpy);
	   SCIPfreeBufferArray(scip, &varbufcpy);
	   SCIPfreeBufferArray(scip, &pscostscores);
	   SCIPfreeBufferArray(scip, &perm);
	   SCIPfreeBufferArray(scip, &randscores);
	   SCIPfreeBufferArray(scip, &redcostscores);
	   SCIPfreeBufferArray(scip, &fixeddistances);
	   SCIPfreeBufferArray(scip, &distances);
	   SCIPfreeBufferArray(scip, &unfixedvars);
	   SCIPfreeBufferArray(scip, &isfixedvar);
	
	   *success = TRUE;
	
	   return SCIP_OKAY;
	}
	
	/** call variable fixing callback for this neighborhood and orchestrate additional variable fixings, if necessary */
	static
	SCIP_RETCODE neighborhoodFixVariables(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   NH*                   neighborhood,       /**< neighborhood data structure */
	   SCIP_VAR**            varbuf,             /**< buffer array to keep variables that should be fixed */
	   SCIP_Real*            valbuf,             /**< buffer array to keep fixing values */
	   int*                  nfixings,           /**< pointer to store the number of variable fixings */
	   SCIP_RESULT*          result              /**< pointer to store the result of the fixing operation */
	   )
	{
	   int ntargetfixings;
	   int nmaxfixings;
	   int nminfixings;
	   int nbinintvars;
	
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	   assert(varbuf != NULL);
	   assert(valbuf != NULL);
	   assert(nfixings != NULL);
	   assert(result != NULL);
	
	   *nfixings = 0;
	
	   *result = SCIP_DIDNOTRUN;
	   ntargetfixings = (int)(neighborhood->fixingrate.targetfixingrate * (SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip)));
	
	   if( neighborhood->varfixings != NULL )
	   {
	      SCIP_CALL( neighborhood->varfixings(scip, neighborhood, varbuf, valbuf, nfixings, result) );
	
	      if( *result != SCIP_SUCCESS )
	         return SCIP_OKAY;
	   }
	   else if( ntargetfixings == 0 )
	   {
	      *result = SCIP_SUCCESS;
	
	      return SCIP_OKAY;
	   }
	
	   /* compute upper and lower target fixing limits using tolerance parameters */
	   assert(neighborhood->varfixings == NULL || *result != SCIP_DIDNOTRUN);
	   nbinintvars = SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip);
	   ntargetfixings = (int)(neighborhood->fixingrate.targetfixingrate * nbinintvars);
	   nminfixings = (int)((neighborhood->fixingrate.targetfixingrate - heurdata->fixtol) * nbinintvars);
	   nminfixings = MAX(nminfixings, 0);
	   nmaxfixings = (int)((neighborhood->fixingrate.targetfixingrate + heurdata->unfixtol) * nbinintvars);
	   nmaxfixings = MIN(nmaxfixings, nbinintvars);
	
	   SCIPdebugMsg(scip, "Neighborhood Fixings/Target: %d / %d <= %d <= %d\n",*nfixings, nminfixings, ntargetfixings, nmaxfixings);
	
	   /* if too few fixings, use a strategy to select more variable fixings: randomized, LP graph, ReducedCost based, mix */
	   if( (*result == SCIP_SUCCESS || *result == SCIP_DIDNOTRUN) && (*nfixings < nminfixings) )
	   {
	      SCIP_Bool success;
	      SCIP_SOL* refsol;
	
	      /* get reference solution from neighborhood */
	      SCIP_CALL( neighborhoodGetRefsol(scip, neighborhood, &refsol) );
	
	      /* try to fix more variables based on the reference solution */
	      if( refsol != NULL )
	      {
	         SCIP_CALL( alnsFixMoreVariables(scip, heurdata, refsol, varbuf, valbuf, nfixings, ntargetfixings, &success) );
	      }
	      else
	         success = FALSE;
	
	      if( success )
	         *result = SCIP_SUCCESS;
	      else if( *result == SCIP_SUCCESS )
	         *result = SCIP_DIDNOTFIND;
	      else
	         *result = SCIP_DIDNOTRUN;
	
	      SCIPdebugMsg(scip, "After additional fixings: %d / %d\n",*nfixings, ntargetfixings);
	   }
	   else if( (SCIP_Real)(*nfixings) > nmaxfixings )
	   {
	      SCIP_Bool success;
	
	      SCIP_CALL( alnsUnfixVariables(scip, heurdata, varbuf, valbuf, nfixings, ntargetfixings, &success) );
	
	      assert(success);
	      *result = SCIP_SUCCESS;
	      SCIPdebugMsg(scip, "Unfixed variables, fixed variables remaining: %d\n", ntargetfixings);
	   }
	   else
	   {
	      SCIPdebugMsg(scip, "No additional fixings performed\n");
	   }
	
	   return SCIP_OKAY;
	}
	
	/** change the sub-SCIP by restricting variable domains, changing objective coefficients, or adding constraints */
	static
	SCIP_RETCODE neighborhoodChangeSubscip(
	   SCIP*                 sourcescip,         /**< source SCIP data structure */
	   SCIP*                 targetscip,         /**< target SCIP data structure */
	   NH*                   neighborhood,       /**< neighborhood */
	   SCIP_VAR**            targetvars,         /**< array of target SCIP variables aligned with source SCIP variables */
	   int*                  ndomchgs,           /**< pointer to store the number of variable domain changes */
	   int*                  nchgobjs,           /**< pointer to store the number of changed objective coefficients */
	   int*                  naddedconss,        /**< pointer to store the number of added constraints */
	   SCIP_Bool*            success             /**< pointer to store whether the sub-SCIP has been successfully modified */
	   )
	{
	   assert(sourcescip != NULL);
	   assert(targetscip != NULL);
	   assert(neighborhood != NULL);
	   assert(targetvars != NULL);
	   assert(ndomchgs != NULL);
	   assert(nchgobjs != NULL);
	   assert(naddedconss != NULL);
	   assert(success != NULL);
	
	   *success = FALSE;
	   *ndomchgs = 0;
	   *nchgobjs = 0;
	   *naddedconss = 0;
	
	   /* call the change sub-SCIP callback of the neighborhood */
	   if( neighborhood->changesubscip != NULL )
	   {
	      SCIP_CALL( neighborhood->changesubscip(sourcescip, targetscip, neighborhood, targetvars, ndomchgs, nchgobjs, naddedconss, success) );
	   }
	   else
	   {
	      *success = TRUE;
	   }
	
	   return SCIP_OKAY;
	}
	
	/** set sub-SCIP solving limits */
	static
	SCIP_RETCODE setLimits(
	   SCIP*                 subscip,            /**< SCIP data structure */
	   SOLVELIMITS*          solvelimits         /**< pointer to solving limits data structure */
	   )
	{
	   assert(subscip != NULL);
	   assert(solvelimits != NULL);
	
	   assert(solvelimits->nodelimit >= solvelimits->stallnodes);
	
	   SCIP_CALL( SCIPsetLongintParam(subscip, "limits/nodes", solvelimits->nodelimit) );
	   SCIP_CALL( SCIPsetLongintParam(subscip, "limits/stallnodes", solvelimits->stallnodes));
	   SCIP_CALL( SCIPsetRealParam(subscip, "limits/time", solvelimits->timelimit) );
	   SCIP_CALL( SCIPsetRealParam(subscip, "limits/memory", solvelimits->memorylimit) );
	
	   return SCIP_OKAY;
	}
	
	/** determine limits for a sub-SCIP */
	static
	SCIP_RETCODE determineLimits(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEUR*            heur,               /**< this heuristic */
	   SOLVELIMITS*          solvelimits,        /**< pointer to solving limits data structure */
	   SCIP_Bool*            runagain            /**< can we solve another sub-SCIP with these limits */
	   )
	{
	   SCIP_HEURDATA* heurdata;
	   SCIP_Real initfactor;
	   SCIP_Real nodesquot;
	   SCIP_Bool avoidmemout;
	
	   assert(scip != NULL);
	   assert(heur != NULL);
	   assert(solvelimits != NULL);
	   assert(runagain != NULL);
	
	   heurdata = SCIPheurGetData(heur);
	
	   /* check whether there is enough time and memory left */
	   SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &solvelimits->timelimit) );
	   if( ! SCIPisInfinity(scip, solvelimits->timelimit) )
	      solvelimits->timelimit -= SCIPgetSolvingTime(scip);
	   SCIP_CALL( SCIPgetRealParam(scip, "limits/memory", &solvelimits->memorylimit) );
	   SCIP_CALL( SCIPgetBoolParam(scip, "misc/avoidmemout", &avoidmemout) );
	
	   /* substract the memory already used by the main SCIP and the estimated memory usage of external software */
	   if( ! SCIPisInfinity(scip, solvelimits->memorylimit) )
	   {
	      solvelimits->memorylimit -= SCIPgetMemUsed(scip)/1048576.0;
	      solvelimits->memorylimit -= SCIPgetMemExternEstim(scip)/1048576.0;
	   }
	
	   /* abort if no time is left or not enough memory (we don't abort in this case if misc_avoidmemout == FALSE)
	   * to create a copy of SCIP, including external memory usage */
	   if( solvelimits->timelimit <= 0.0 || (avoidmemout && solvelimits->memorylimit <= 2.0*SCIPgetMemExternEstim(scip)/1048576.0) )
	      *runagain = FALSE;
	
	   nodesquot = heurdata->nodesquot;
	
	   /* if the heuristic is used to measure all rewards, it will always be penalized here */
	   if( heurdata->rewardfile == NULL )
	      nodesquot *= (SCIPheurGetNBestSolsFound(heur) + 1.0)/(SCIPheurGetNCalls(heur) + 1.0);
	
	   nodesquot = MAX(nodesquot, heurdata->nodesquotmin);
	
	   /* calculate the search node limit of the heuristic  */
	   solvelimits->stallnodes = (SCIP_Longint)(nodesquot * SCIPgetNNodes(scip));
	   solvelimits->stallnodes += heurdata->nodesoffset;
	   solvelimits->stallnodes -= heurdata->usednodes;
	   solvelimits->stallnodes -= 100 * SCIPheurGetNCalls(heur);
	   solvelimits->stallnodes = MIN(heurdata->maxnodes, solvelimits->stallnodes);
	
	   /* use a smaller budget if not all neighborhoods have been initialized yet */
	   assert(heurdata->ninitneighborhoods >= 0);
	   initfactor = (heurdata->nactiveneighborhoods - heurdata->ninitneighborhoods + 1.0) / (heurdata->nactiveneighborhoods + 1.0);
	   solvelimits->stallnodes = (SCIP_Longint)(solvelimits->stallnodes * initfactor);
	   solvelimits->nodelimit = (SCIP_Longint)(heurdata->maxnodes);
	
	   /* check whether we have enough nodes left to call subproblem solving */
	   if( solvelimits->stallnodes < heurdata->targetnodes )
	      *runagain = FALSE;
	
	   return SCIP_OKAY;
	}
	
	/** return the bandit algorithm that should be used */
	static
	SCIP_BANDIT* getBandit(
	   SCIP_HEURDATA*        heurdata            /**< heuristic data of the ALNS neighborhood */
	   )
	{
	   assert(heurdata != NULL);
	   return heurdata->bandit;
	}
	
	/** select a neighborhood depending on the selected bandit algorithm */
	static
	SCIP_RETCODE selectNeighborhood(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   int*                  neighborhoodidx     /**< pointer to store the selected neighborhood index */
	   )
	{
	   SCIP_BANDIT* bandit;
	   assert(scip != NULL);
	   assert(heurdata != NULL);
	   assert(neighborhoodidx != NULL);
	
	   *neighborhoodidx = -1;
	
	   bandit = getBandit(heurdata);
	
	   SCIP_CALL( SCIPbanditSelect(bandit, neighborhoodidx) );
	   assert(*neighborhoodidx >= 0);
	
	   return SCIP_OKAY;
	}
	
	/** Calculate reward based on the selected reward measure */
	static
	SCIP_RETCODE getReward(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   NH_STATS*             runstats,           /**< run statistics */
	   SCIP_Real*            rewardptr           /**< array to store the computed rewards, total and individual */
	   )
	{
	   SCIP_Real reward = 0.0;
	   SCIP_Real effort;
	   int ndiscretevars;
	
	   memset(rewardptr, 0, sizeof(*rewardptr)*(int)NREWARDTYPES);
	
	   assert(rewardptr != NULL);
	   assert(runstats->usednodes >= 0);
	   assert(runstats->nfixings >= 0);
	
	   effort = runstats->usednodes / 100.0;
	
	   ndiscretevars = SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip);
	   /* assume that every fixed variable linearly reduces the subproblem complexity */
	   if( ndiscretevars > 0 )
	   {
	      effort = (1.0 - (runstats->nfixings / (SCIP_Real)ndiscretevars)) * effort;
	   }
	   assert(rewardptr != NULL);
	
	   /* a positive reward is only assigned if a new incumbent solution was found */
	   if( runstats->nbestsolsfound > 0 )
	   {
	      SCIP_Real rewardcontrol = heurdata->rewardcontrol;
	
	      SCIP_Real lb;
	      SCIP_Real ub;
	
	      /* the indicator function is simply 1.0 */
	      rewardptr[REWARDTYPE_BESTSOL] = 1.0;
	      rewardptr[REWARDTYPE_NOSOLPENALTY] = 1.0;
	
	      ub = runstats->newupperbound;
	      lb = SCIPgetLowerbound(scip);
	
	      /* compute the closed gap reward */
	      if( SCIPisEQ(scip, ub, lb) || SCIPisInfinity(scip, runstats->oldupperbound) )
	         rewardptr[REWARDTYPE_CLOSEDGAP] = 1.0;
	      else
	      {
	         rewardptr[REWARDTYPE_CLOSEDGAP] = (runstats->oldupperbound - ub) / (runstats->oldupperbound - lb);
	      }
	
	      /* the reward is a convex combination of the best solution reward and the reward for the closed gap */
	      reward = rewardcontrol * rewardptr[REWARDTYPE_BESTSOL] + (1.0 - rewardcontrol) * rewardptr[REWARDTYPE_CLOSEDGAP];
	
	      /* optionally, scale the reward by the involved effort */
	      if( heurdata->scalebyeffort )
	         reward /= (effort + 1.0);
	
	      /* add the baseline and rescale the reward into the interval [baseline, 1.0] */
	      reward = heurdata->rewardbaseline + (1.0 - heurdata->rewardbaseline) * reward;
	   }
	   else
	   {
	      /* linearly decrease the reward based on the number of nodes spent */
	      SCIP_Real maxeffort = heurdata->targetnodes;
	      SCIP_Real usednodes = runstats->usednodes;
	
	      if( ndiscretevars > 0 )
	         usednodes *= (1.0 - (runstats->nfixings / (SCIP_Real)ndiscretevars));
	
	      rewardptr[REWARDTYPE_NOSOLPENALTY] = 1 - (usednodes / maxeffort);
	      rewardptr[REWARDTYPE_NOSOLPENALTY] = MAX(0.0, rewardptr[REWARDTYPE_NOSOLPENALTY]);
	      reward = heurdata->rewardbaseline * rewardptr[REWARDTYPE_NOSOLPENALTY];
	   }
	
	   rewardptr[REWARDTYPE_TOTAL] = reward;
	
	   return SCIP_OKAY;
	}
	
	/** update internal bandit algorithm statistics for future draws */
	static
	SCIP_RETCODE updateBanditAlgorithm(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata,           /**< heuristic data of the ALNS neighborhood */
	   SCIP_Real             reward,             /**< measured reward */
	   int                   neighborhoodidx     /**< the neighborhood that was chosen */
	   )
	{
	   SCIP_BANDIT* bandit;
	   assert(scip != NULL);
	   assert(heurdata != NULL);
	   assert(neighborhoodidx >= 0);
	   assert(neighborhoodidx < heurdata->nactiveneighborhoods);
	
	   bandit = getBandit(heurdata);
	
	   SCIPdebugMsg(scip, "Rewarding bandit algorithm action %d with reward %.2f\n", neighborhoodidx, reward);
	   SCIP_CALL( SCIPbanditUpdate(bandit, neighborhoodidx, reward) );
	
	   return SCIP_OKAY;
	}
	
	/** set up the sub-SCIP parameters, objective cutoff, and solution limits */
	static
	SCIP_RETCODE setupSubScip(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP*                 subscip,            /**< sub-SCIP data structure */
	   SCIP_VAR**            subvars,            /**< array of sub-SCIP variables in the order of the main SCIP */
	   SOLVELIMITS*          solvelimits,        /**< pointer to solving limits data structure */
	   SCIP_HEUR*            heur,               /**< this heuristic */
	   SCIP_Bool             objchgd             /**< did the objective change between the source and the target SCIP? */
	   )
	{
	   SCIP_HEURDATA* heurdata;
	   SCIP_Real cutoff;
	   SCIP_Real upperbound;
	
	   heurdata = SCIPheurGetData(heur);
	
	   /* do not abort subproblem on CTRL-C */
	   SCIP_CALL( SCIPsetBoolParam(subscip, "misc/catchctrlc", FALSE) );
	
	   /* disable output to console unless we are in debug mode */
	   SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 0) );
	
	   /* disable statistic timing inside sub SCIP */
	   SCIP_CALL( SCIPsetBoolParam(subscip, "timing/statistictiming", FALSE) );
	
	#ifdef ALNS_SUBSCIPOUTPUT
	   SCIP_CALL( SCIPsetIntParam(subscip, "display/verblevel", 5) );
	   SCIP_CALL( SCIPsetIntParam(subscip, "display/freq", 1) );
	   /* enable statistic timing inside sub SCIP */
	      SCIP_CALL( SCIPsetBoolParam(subscip, "timing/statistictiming", TRUE) );
	#endif
	
	   SCIP_CALL( SCIPsetIntParam(subscip, "limits/bestsol", heurdata->nsolslim) );
	
	   /* forbid recursive call of heuristics and separators solving subMIPs */
	   if( ! heurdata->usesubscipheurs )
	   {
	      SCIP_CALL( SCIPsetSubscipsOff(subscip, TRUE) );
	   }
	
	   /* disable cutting plane separation */
	   SCIP_CALL( SCIPsetSeparating(subscip, SCIP_PARAMSETTING_OFF, TRUE) );
	
	   /* disable expensive presolving */
	   SCIP_CALL( SCIPsetPresolving(subscip, SCIP_PARAMSETTING_FAST, TRUE) );
	
	   /* use best estimate node selection */
	   if( SCIPfindNodesel(subscip, "estimate") != NULL && ! SCIPisParamFixed(subscip, "nodeselection/estimate/stdpriority") )
	   {
	      SCIP_CALL( SCIPsetIntParam(subscip, "nodeselection/estimate/stdpriority", INT_MAX/4) );
	   }
	
	   /* use inference branching */

	   if( SCIPfindBranchrule(subscip, "inference") != NULL && ! SCIPisParamFixed(subscip, "branching/inference/priority") )
	   {
	      SCIP_CALL( SCIPsetIntParam(subscip, "branching/inference/priority", INT_MAX/4) );
	   }
	
	   /* enable conflict analysis and restrict conflict pool */
	   if( ! SCIPisParamFixed(subscip, "conflict/enable") )
	   {
	      SCIP_CALL( SCIPsetBoolParam(subscip, "conflict/enable", TRUE) );
	   }
	
	   if( !SCIPisParamFixed(subscip, "conflict/useboundlp") )
	   {
	      SCIP_CALL( SCIPsetCharParam(subscip, "conflict/useboundlp", 'o') );
	   }
	
	   if( ! SCIPisParamFixed(subscip, "conflict/maxstoresize") )
	   {
	      SCIP_CALL( SCIPsetIntParam(subscip, "conflict/maxstoresize", 100) );
	   }
	
	   /* speed up sub-SCIP by not checking dual LP feasibility */
	   SCIP_CALL( SCIPsetBoolParam(subscip, "lp/checkdualfeas", FALSE) );
	
	   /* add an objective cutoff */
	   if( ! SCIPisInfinity(scip, SCIPgetUpperbound(scip)) )
	   {
	      upperbound = SCIPgetUpperbound(scip) - SCIPsumepsilon(scip);
	      if( ! SCIPisInfinity(scip, -1.0 * SCIPgetLowerbound(scip)) )
	      {
	         cutoff = (1 - heurdata->minimprove) * SCIPgetUpperbound(scip)
	                            + heurdata->minimprove * SCIPgetLowerbound(scip);
	      }
	      else
	      {
	         if( SCIPgetUpperbound(scip) >= 0 )
	            cutoff = (1 - heurdata->minimprove) * SCIPgetUpperbound(scip);
	         else
	            cutoff = (1 + heurdata->minimprove) * SCIPgetUpperbound(scip);
	      }
	      cutoff = MIN(upperbound, cutoff);
	
	      if( SCIPisObjIntegral(scip) )
	         cutoff = SCIPfloor(scip, cutoff);
	
	      SCIPdebugMsg(scip, "Sub-SCIP cutoff: %15.9" SCIP_REAL_FORMAT " (%15.9" SCIP_REAL_FORMAT " in original space)\n",
	         cutoff, SCIPretransformObj(scip, cutoff));
	
	      /* if the objective changed between the source and the target SCIP, encode the cutoff as a constraint */
	      if( ! objchgd )
	      {
	         SCIP_CALL(SCIPsetObjlimit(subscip, cutoff));
	
	         SCIPdebugMsg(scip, "Cutoff added as Objective Limit\n");
	      }
	      else
	      {
	         SCIP_CONS* objcons;
	         int nvars;
	         SCIP_VAR** vars;
	         int i;
	
	         vars = SCIPgetVars(scip);
	         nvars = SCIPgetNVars(scip);
	
	         SCIP_CALL( SCIPcreateConsLinear(subscip, &objcons, "objbound_of_origscip", 0, NULL, NULL, -SCIPinfinity(subscip), cutoff,
	                     TRUE, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE) );
	         for( i = 0; i < nvars; ++i)
	         {
	            if( ! SCIPisFeasZero(subscip, SCIPvarGetObj(vars[i])) )
	            {
	               assert(subvars[i] != NULL);
	               SCIP_CALL( SCIPaddCoefLinear(subscip, objcons, subvars[i], SCIPvarGetObj(vars[i])) );
	            }
	         }
	         SCIP_CALL( SCIPaddCons(subscip, objcons) );
	         SCIP_CALL( SCIPreleaseCons(subscip, &objcons) );
	
	         SCIPdebugMsg(scip, "Cutoff added as constraint\n");
	      }
	   }
	
	   /* set solve limits for sub-SCIP */
	   SCIP_CALL( setLimits(subscip, solvelimits) );
	
	   /* change random seed of sub-SCIP */
	   if( heurdata->subsciprandseeds )
	   {
	      SCIP_CALL( SCIPsetIntParam(subscip, "randomization/randomseedshift", (int)SCIPheurGetNCalls(heur)) );
	   }
	
	   SCIPdebugMsg(scip, "Solve Limits: %lld (%lld) nodes (stall nodes), %.1f sec., %d sols\n",
	         solvelimits->nodelimit, solvelimits->stallnodes, solvelimits->timelimit, heurdata->nsolslim);
	
	   return SCIP_OKAY;
	}
	
	/** execution method of primal heuristic */
	static
	SCIP_DECL_HEUREXEC(heurExecAlns)
	{  /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	   SCIP_VAR** varbuf;
	   SCIP_Real* valbuf;
	   SCIP_VAR** vars;
	   SCIP_VAR** subvars;
	   NH_STATS runstats[NNEIGHBORHOODS];
	   SCIP_STATUS subscipstatus[NNEIGHBORHOODS];
	   SCIP* subscip = NULL;
	
	   int nfixings;
	   int nvars;
	   int neighborhoodidx;
	   int ntries;
	   SCIP_Bool tryagain;
	   NH* neighborhood;
	   SOLVELIMITS solvelimits;
	   SCIP_Bool success;
	   SCIP_Bool run;
	   SCIP_Bool allrewardsmode;
	   SCIP_Real rewards[NNEIGHBORHOODS][NREWARDTYPES] = {{0}};
	   int banditidx;
	
	   int i;
	
	   heurdata = SCIPheurGetData(heur);
	   assert(heurdata != NULL);
	
	   *result = SCIP_DIDNOTRUN;
	
	   if( heurdata->nactiveneighborhoods == 0 )
	      return SCIP_OKAY;
	
	   /* we only allow to run multiple times at a node during the root */
	   if( (heurtiming & SCIP_HEURTIMING_DURINGLPLOOP) && (SCIPgetDepth(scip) > 0 || !heurdata->initduringroot) )
	      return SCIP_OKAY;
	
	   /* update internal incumbent solution */
	   if( SCIPgetBestSol(scip) != heurdata->lastcallsol )
	   {
	      heurdata->lastcallsol = SCIPgetBestSol(scip);
	      heurdata->firstcallthissol = SCIPheurGetNCalls(heur);
	   }
	
	   /* do not run more than a user-defined number of times on each incumbent (-1: no limit) */
	   if( heurdata->maxcallssamesol != -1 )
	   {
	      SCIP_Longint samesollimit = (heurdata->maxcallssamesol > 0) ?
	         heurdata->maxcallssamesol :
	         heurdata->nactiveneighborhoods;
	
	      if( SCIPheurGetNCalls(heur) - heurdata->firstcallthissol >= samesollimit )
	      {
	         SCIPdebugMsg(scip, "Heuristic already called %" SCIP_LONGINT_FORMAT " times on current incumbent\n", SCIPheurGetNCalls(heur) - heurdata->firstcallthissol);
	         return SCIP_OKAY;
	      }
	   }
	
	   /* wait for a sufficient number of nodes since last incumbent solution */
	   if( SCIPgetDepth(scip) > 0 && SCIPgetBestSol(scip) != NULL
	      && (SCIPgetNNodes(scip) - SCIPsolGetNodenum(SCIPgetBestSol(scip))) < heurdata->waitingnodes )
	   {
	      SCIPdebugMsg(scip, "Waiting nodes not satisfied\n");
	      return SCIP_OKAY;
	   }
	
	   run = TRUE;
	   /* check if budget allows a run of the next selected neighborhood */
	   SCIP_CALL( determineLimits(scip, heur, &solvelimits, &run) );
	   SCIPdebugMsg(scip, "Budget check: %" SCIP_LONGINT_FORMAT " (%" SCIP_LONGINT_FORMAT ") %s\n", solvelimits.nodelimit, heurdata->targetnodes, run ? "passed" : "must wait");
	
	   if( ! run )
	      return SCIP_OKAY;
	
	   /* delay the heuristic if local reduced costs should be used for generic variable unfixing */
	   if( heurdata->uselocalredcost && (nodeinfeasible || ! SCIPhasCurrentNodeLP(scip) || SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL) )
	   {
	      *result = SCIP_DELAYED;
	
	      return SCIP_OKAY;
	   }
	
	   allrewardsmode = heurdata->rewardfile != NULL;
	
	   /* apply some other rules for a fair all rewards mode; in normal execution mode, neighborhoods are iterated through */
	   if( allrewardsmode )
	   {
	      /* most neighborhoods require an incumbent solution */
	      if( SCIPgetNSols(scip) < 2 )
	      {
	         SCIPdebugMsg(scip, "Not enough solutions for all rewards mode\n");
	         return SCIP_OKAY;
	      }
	
	      /* if the node is infeasible, or has no LP solution, which is required by some neighborhoods
	       * if we are not in all rewards mode, the neighborhoods delay themselves individually
	       */
	      if( nodeinfeasible || ! SCIPhasCurrentNodeLP(scip) || SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
	      {
	         SCIPdebugMsg(scip, "Delay ALNS heuristic until a feasible node with optimally solved LP relaxation\n");
	         *result = SCIP_DELAYED;
	         return SCIP_OKAY;
	      }
	   }
	
	   /* use the neighborhood that requested a delay or select the next neighborhood to run based on the selected bandit algorithm */
	   if( heurdata->currneighborhood >= 0 )
	   {
	      assert(! allrewardsmode);
	      banditidx = heurdata->currneighborhood;
	      SCIPdebugMsg(scip, "Select delayed neighborhood %d (was delayed %d times)\n", banditidx, heurdata->ndelayedcalls);
	   }
	   else
	   {
	      SCIP_CALL( selectNeighborhood(scip, heurdata, &banditidx) );
	      SCIPdebugMsg(scip, "Selected neighborhood %d with bandit algorithm\n", banditidx);
	   }
	
	   /* in all rewards mode, we simply loop over all heuristics */
	   if( ! allrewardsmode )
	      neighborhoodidx = banditidx;
	   else
	      neighborhoodidx = 0;
	
	   assert(0 <= neighborhoodidx && neighborhoodidx < NNEIGHBORHOODS);
	   assert(heurdata->nactiveneighborhoods > neighborhoodidx);
	
	   /* allocate memory for variable fixings buffer */
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, NULL, NULL, NULL, NULL) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &varbuf, nvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &valbuf, nvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &subvars, nvars) );
	
	   /* initialize neighborhood statistics for a run */
	   ntries = 1;
	   do
	   {
	      SCIP_HASHMAP* varmapf;
	      SCIP_EVENTHDLR* eventhdlr;
	      SCIP_EVENTDATA eventdata;
	      char probnamesuffix[SCIP_MAXSTRLEN];
	      SCIP_Real allfixingrate;
	      int ndomchgs;
	      int nchgobjs;
	      int naddedconss;
	      int v;
	      SCIP_RETCODE retcode;
	      SCIP_RESULT fixresult;
	
	      tryagain = FALSE;
	      neighborhood = heurdata->neighborhoods[neighborhoodidx];
	      SCIPdebugMsg(scip, "Running '%s' neighborhood %d\n", neighborhood->name, neighborhoodidx);
	
	      initRunStats(scip, &runstats[neighborhoodidx]);
	      rewards[neighborhoodidx][REWARDTYPE_TOTAL] = 0.0;
	
	      subscipstatus[neighborhoodidx] = SCIP_STATUS_UNKNOWN;
	      SCIP_CALL( SCIPstartClock(scip, neighborhood->stats.setupclock) );
	
	      /* determine variable fixings and objective coefficients of this neighborhood */
	      SCIP_CALL( neighborhoodFixVariables(scip, heurdata, neighborhood, varbuf, valbuf, &nfixings, &fixresult) );
	
	      SCIPdebugMsg(scip, "Fix %d/%d variables, result code %d\n", nfixings, nvars,fixresult);
	
	      /* Fixing was not successful, either because the fixing rate was not reached (and no additional variable
	       * prioritization was used), or the neighborhood requested a delay, e.g., because no LP relaxation solution exists
	       * at the current node
	       *
	       * The ALNS heuristic keeps a delayed neighborhood active and delays itself.
	       */
	      if( fixresult != SCIP_SUCCESS )
	      {
	         SCIP_CALL( SCIPstopClock(scip, neighborhood->stats.setupclock) );
	
	         /* to determine all rewards, we cannot delay neighborhoods */
	         if( allrewardsmode )
	         {
	            if( ntries == heurdata->nactiveneighborhoods )
	               break;
	
	            neighborhoodidx = (neighborhoodidx + 1) % heurdata->nactiveneighborhoods;
	            ntries++;
	            tryagain = TRUE;
	
	            continue;
	         }
	
	         /* delay the heuristic along with the selected neighborhood
	          *
	          * if the neighborhood has been delayed for too many consecutive calls, the delay is treated as a failure */
	         if( fixresult == SCIP_DELAYED )
	         {
	            if( heurdata->ndelayedcalls > (SCIPheurGetFreq(heur) / 4 + 1) )
	            {
	               resetCurrentNeighborhood(heurdata);
	
	               /* use SCIP_DIDNOTFIND to penalize the neighborhood with a bad reward */
	               fixresult = SCIP_DIDNOTFIND;
	            }
	            else if( heurdata->currneighborhood == -1 )
	            {
	               heurdata->currneighborhood = neighborhoodidx;
	               heurdata->ndelayedcalls = 1;
	            }
	            else
	            {
	               heurdata->ndelayedcalls++;
	            }
	         }
	
	         if( fixresult == SCIP_DIDNOTRUN )
	         {
	            if( ntries < heurdata->nactiveneighborhoods )
	            {
	               SCIP_CALL( updateBanditAlgorithm(scip, heurdata, 0.0, neighborhoodidx) );
	               SCIP_CALL( selectNeighborhood(scip, heurdata, &neighborhoodidx) );
	               ntries++;
	               tryagain = TRUE;
	
	               SCIPdebugMsg(scip, "Neighborhood cannot run -> try next neighborhood %d\n", neighborhoodidx);
	               continue;
	            }
	            else
	               break;
	         }
	
	         assert(fixresult == SCIP_DIDNOTFIND || fixresult == SCIP_DELAYED);
	         *result = fixresult;
	         break;
	      }
	
	      *result = SCIP_DIDNOTFIND;
	
	      neighborhood->stats.nfixings += nfixings;
	      runstats[neighborhoodidx].nfixings = nfixings;
	
	      SCIP_CALL( SCIPcreate(&subscip) );
	      SCIP_CALL( SCIPhashmapCreate(&varmapf, SCIPblkmem(scip), nvars) );
	      (void) SCIPsnprintf(probnamesuffix, SCIP_MAXSTRLEN, "alns_%s", neighborhood->name);
	
	      /* todo later: run global propagation for this set of fixings */
	      SCIP_CALL( SCIPcopyLargeNeighborhoodSearch(scip, subscip, varmapf, probnamesuffix, varbuf, valbuf, nfixings, FALSE, heurdata->copycuts, &success, NULL) );
	
	      /* store sub-SCIP variables in array for faster access */
	      for( v = 0; v < nvars; ++v )
	      {
	         subvars[v] = (SCIP_VAR*)SCIPhashmapGetImage(varmapf, (void *)vars[v]);
	      }
	
	      SCIPhashmapFree(&varmapf);
	
	      /* let the neighborhood add additional constraints, or restrict domains */
	      SCIP_CALL( neighborhoodChangeSubscip(scip, subscip, neighborhood, subvars, &ndomchgs, &nchgobjs, &naddedconss, &success) );
	
	      if( ! success )
	      {
	         SCIP_CALL( SCIPstopClock(scip, neighborhood->stats.setupclock) );
	
	         if( ! allrewardsmode || ntries == heurdata->nactiveneighborhoods )
	            break;
	
	         neighborhoodidx = (neighborhoodidx + 1) % heurdata->nactiveneighborhoods;
	         ntries++;
	         tryagain = TRUE;
	
	         SCIP_CALL( SCIPfree(&subscip) );
	
	         continue;
	      }
	
	      /* set up sub-SCIP parameters */
	      SCIP_CALL( setupSubScip(scip, subscip, subvars, &solvelimits, heur, nchgobjs > 0) );
	
	      /* copy the necessary data into the event data to create new solutions */
	      eventdata.nodelimit = solvelimits.nodelimit;  /*lint !e644*/
	      eventdata.lplimfac = heurdata->lplimfac;
	      eventdata.heur = heur;
	      eventdata.sourcescip = scip;
	      eventdata.subvars = subvars;
	      eventdata.runstats = &runstats[neighborhoodidx];
	      eventdata.allrewardsmode = allrewardsmode;
	
	      /* include an event handler to transfer solutions into the main SCIP */
	      SCIP_CALL( SCIPincludeEventhdlrBasic(subscip, &eventhdlr, EVENTHDLR_NAME, EVENTHDLR_DESC, eventExecAlns, NULL) );
	
	      /* transform the problem before catching the events */
	      SCIP_CALL( SCIPtransformProb(subscip) );
	      SCIP_CALL( SCIPcatchEvent(subscip, SCIP_EVENTTYPE_ALNS, eventhdlr, &eventdata, NULL) );
	
	      SCIP_CALL( SCIPstopClock(scip, neighborhood->stats.setupclock) );
	
	      SCIP_CALL( SCIPstartClock(scip, neighborhood->stats.submipclock) );
	
	      /* set up sub-SCIP and run presolving */
	      retcode = SCIPpresolve(subscip);
	      if( retcode != SCIP_OKAY )
	      {
	         SCIPwarningMessage(scip, "Error while presolving subproblem in ALNS heuristic; sub-SCIP terminated with code <%d>\n", retcode);
	         SCIP_CALL( SCIPstopClock(scip, neighborhood->stats.submipclock) );
	
	         SCIPABORT();  /*lint --e{527}*/
	         break;
	      }
	
	      /* was presolving successful enough regarding fixings? otherwise, terminate */
	      allfixingrate = (SCIPgetNOrigVars(subscip) - SCIPgetNVars(subscip)) / (SCIP_Real)SCIPgetNOrigVars(subscip);
	
	      /* additional variables added in presolving may lead to the subSCIP having more variables than the original */
	      allfixingrate = MAX(allfixingrate, 0.0);
	
	      if( allfixingrate >= neighborhood->fixingrate.targetfixingrate / 2.0 )
	      {
	         /* run sub-SCIP for the given budget, and collect statistics */
	         SCIP_CALL_ABORT( SCIPsolve(subscip) );
	      }
	      else
	      {
	         SCIPdebugMsg(scip, "Fixed only %.3f of all variables after presolving -> do not solve sub-SCIP\n", allfixingrate);
	      }
	
	#ifdef ALNS_SUBSCIPOUTPUT
	      SCIP_CALL( SCIPprintStatistics(subscip, NULL) );
	#endif
	
	      SCIP_CALL( SCIPstopClock(scip, neighborhood->stats.submipclock) );
	
	      /* update statistics based on the sub-SCIP run results */
	      updateRunStats(&runstats[neighborhoodidx], subscip);
	      subscipstatus[neighborhoodidx] = SCIPgetStatus(subscip);
	      SCIPdebugMsg(scip, "Status of sub-SCIP run: %d\n", subscipstatus[neighborhoodidx]);
	
	      SCIP_CALL( getReward(scip, heurdata, &runstats[neighborhoodidx], rewards[neighborhoodidx]) );
	
	      /* in all rewards mode, continue with the next neighborhood */
	      if( allrewardsmode && ntries < heurdata->nactiveneighborhoods )
	      {
	         neighborhoodidx = (neighborhoodidx + 1) % heurdata->nactiveneighborhoods;
	         ntries++;
	         tryagain = TRUE;
	
	         SCIP_CALL( SCIPfree(&subscip) );
	      }
	   }
	   while( tryagain && ! SCIPisStopped(scip) );
	
	   if( subscip != NULL )
	   {
	      SCIP_CALL( SCIPfree(&subscip) );
	   }
	
	   SCIPfreeBufferArray(scip, &subvars);
	   SCIPfreeBufferArray(scip, &valbuf);
	   SCIPfreeBufferArray(scip, &varbuf);
	
	   /* update bandit index that may have changed unless we are in all rewards mode */
	   if( ! allrewardsmode )
	      banditidx = neighborhoodidx;
	
	   if( *result != SCIP_DELAYED )
	   {
	      /* decrease the number of neighborhoods that have not been initialized */
	      if( neighborhood->stats.nruns == 0 )
	         --heurdata->ninitneighborhoods;
	
	      heurdata->usednodes += runstats[banditidx].usednodes;
	
	      /* determine the success of this neighborhood, and update the target fixing rate for the next time */
	      updateNeighborhoodStats(&runstats[banditidx], heurdata->neighborhoods[banditidx], subscipstatus[banditidx]);
	
	      /* adjust the fixing rate for this neighborhood
	       * make no adjustments in all rewards mode, because this only affects 1 of 8 heuristics
	       */
	      if( heurdata->adjustfixingrate && ! allrewardsmode )
	      {
	         SCIPdebugMsg(scip, "Update fixing rate: %.2f\n", heurdata->neighborhoods[banditidx]->fixingrate.targetfixingrate);
	         updateFixingRate(heurdata->neighborhoods[banditidx], subscipstatus[banditidx], &runstats[banditidx]);
	         SCIPdebugMsg(scip, "New fixing rate: %.2f\n", heurdata->neighborhoods[banditidx]->fixingrate.targetfixingrate);
	      }
	      /* similarly, update the minimum improvement for the ALNS heuristic */
	      if( heurdata->adjustminimprove )
	      {
	         SCIPdebugMsg(scip, "Update Minimum Improvement: %.4f\n", heurdata->minimprove);
	         updateMinimumImprovement(heurdata, subscipstatus[banditidx], &runstats[banditidx]);
	         SCIPdebugMsg(scip, "--> %.4f\n", heurdata->minimprove);
	      }
	
	      /* update the target node limit based on the status of the selected algorithm */
	      if( heurdata->adjusttargetnodes && SCIPheurGetNCalls(heur) >= heurdata->nactiveneighborhoods )
	      {
	         updateTargetNodeLimit(heurdata, &runstats[banditidx], subscipstatus[banditidx]);
	      }
	
	      /* update the bandit algorithm by the measured reward */
	      SCIP_CALL( updateBanditAlgorithm(scip, heurdata, rewards[banditidx][REWARDTYPE_TOTAL], banditidx) );
	
	      resetCurrentNeighborhood(heurdata);
	   }
	
	   /* write single, measured rewards and the bandit index to the reward file */
	   if( allrewardsmode )
	   {
	      int j;
	      for( j = 0; j < (int)NREWARDTYPES; j++ )
	         for( i = 0; i < heurdata->nactiveneighborhoods; ++i )
	            fprintf(heurdata->rewardfile, "%.4f,", rewards[i][j]);
	
	      fprintf(heurdata->rewardfile, "%d\n", banditidx);
	   }
	
	   return SCIP_OKAY;
	}
	
	/** callback to collect variable fixings of RENS */
	static
	DECL_VARFIXINGS(varFixingsRens)
	{  /*lint --e{715}*/
	   int nbinvars;
	   int nintvars;
	   SCIP_VAR** vars;
	   int i;
	   int *fracidx = NULL;
	   SCIP_Real* frac = NULL;
	   int nfracs;
	
	   assert(scip != NULL);
	   assert(varbuf != NULL);
	   assert(nfixings != NULL);
	   assert(valbuf != NULL);
	
	   *result = SCIP_DELAYED;
	
	   if( ! SCIPhasCurrentNodeLP(scip) )
	      return SCIP_OKAY;
	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
	      return SCIP_OKAY;
	
	   *result = SCIP_DIDNOTRUN;
	
	   /* get variable information */
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	
	   /* return if no binary or integer variables are present */
	   if( nbinvars + nintvars == 0 )
	      return SCIP_OKAY;
	
	   SCIP_CALL( SCIPallocBufferArray(scip, &fracidx, nbinvars + nintvars) );
	   SCIP_CALL( SCIPallocBufferArray(scip, &frac, nbinvars + nintvars) );
	
	   /* loop over binary and integer variables; determine those that should be fixed in the sub-SCIP */
	   for( nfracs = 0, i = 0; i < nbinvars + nintvars; ++i )
	   {
	      SCIP_VAR* var = vars[i];
	      SCIP_Real lpsolval = SCIPvarGetLPSol(var);
	      assert((i < nbinvars && SCIPvarIsBinary(var)) || (i >= nbinvars && SCIPvarIsIntegral(var)));
	
	      /* fix all binary and integer variables with integer LP solution value */
	      if( SCIPisFeasIntegral(scip, lpsolval) )
	      {
	         tryAdd2variableBuffer(scip, var, lpsolval, varbuf, valbuf, nfixings, TRUE);
	      }
	      else
	      {
	         frac[nfracs] = SCIPfrac(scip, lpsolval);
	         frac[nfracs] = MIN(frac[nfracs], 1.0 - frac[nfracs]);
	         fracidx[nfracs++] = i;
	      }
	   }
	
	   /* do some additional fixing */
	   if( *nfixings < neighborhood->fixingrate.targetfixingrate * (nbinvars + nintvars) && nfracs > 0 )
	   {
	      SCIPsortDownRealInt(frac, fracidx, nfracs);
	
	      /* prefer variables that are almost integer */
	      for( i = 0; i < nfracs && *nfixings < neighborhood->fixingrate.targetfixingrate * (nbinvars + nintvars); i++ )
	      {
	         tryAdd2variableBuffer(scip, vars[fracidx[i]], SCIPround(scip, SCIPvarGetLPSol(vars[fracidx[i]])), varbuf, valbuf, nfixings, TRUE);
	      }
	   }
	
	   SCIPfreeBufferArray(scip, &frac);
	   SCIPfreeBufferArray(scip, &fracidx);
	
	   *result = SCIP_SUCCESS;
	
	   return SCIP_OKAY;
	}
	
	/** callback for RENS subproblem changes */
	static
	DECL_CHANGESUBSCIP(changeSubscipRens)
	{  /*lint --e{715}*/
	   SCIP_VAR** vars;
	   int nintvars;
	   int nbinvars;
	   int i;
	
	   assert(SCIPhasCurrentNodeLP(sourcescip));
	   assert(SCIPgetLPSolstat(sourcescip) == SCIP_LPSOLSTAT_OPTIMAL);
	
	   /* get variable information */
	   SCIP_CALL( SCIPgetVarsData(sourcescip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	
	   /* restrict bounds of integer variables with fractional solution value */
	   for( i = nbinvars; i < nbinvars + nintvars; ++i )
	   {
	      SCIP_VAR* var = vars[i];
	      SCIP_Real lpsolval = SCIPgetSolVal(sourcescip, NULL, var);
	
	      if( subvars[i] == NULL )
	         continue;
	
	      if( ! SCIPisFeasIntegral(sourcescip, lpsolval) )
	      {
	         SCIP_Real newlb = SCIPfloor(sourcescip, lpsolval);
	         SCIP_Real newub = newlb + 1.0;
	
	         /* only count this as a domain change if the new lower and upper bound are a further restriction */
	         if( newlb > SCIPvarGetLbGlobal(subvars[i]) + 0.5 || newub < SCIPvarGetUbGlobal(subvars[i]) - 0.5 )
	         {
	            SCIP_CALL( SCIPchgVarLbGlobal(targetscip, subvars[i], newlb) );
	            SCIP_CALL( SCIPchgVarUbGlobal(targetscip, subvars[i], newub) );
	            (*ndomchgs)++;
	         }
	      }
	   }
	
	   *success = TRUE;
	
	   return SCIP_OKAY;
	}
	
	/** collect fixings by matching solution values in a collection of solutions for all binary and integer variables,
	 *  or for a custom set of variables
	 */
	static
	SCIP_RETCODE fixMatchingSolutionValues(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_SOL**            sols,               /**< array of 2 or more solutions. It is okay for the array to contain one element
	                                              *   equal to NULL to represent the current LP solution */
	   int                   nsols,              /**< number of solutions in the array */
	   SCIP_VAR**            vars,               /**< variable array for which solution values must agree */
	   int                   nvars,              /**< number of variables, or -1 for all binary and integer variables */
	   SCIP_VAR**            varbuf,             /**< buffer storage for variable fixings */
	   SCIP_Real*            valbuf,             /**< buffer storage for fixing values */
	   int*                  nfixings            /**< pointer to store the number of fixings */
	   )
	{
	   int v;
	   int nbinintvars;
	   SCIP_SOL* firstsol;
	
	   assert(scip != NULL);
	   assert(sols != NULL);
	   assert(nsols >= 2);
	   assert(varbuf != NULL);
	   assert(valbuf != NULL);
	   assert(nfixings != NULL);
	   assert(*nfixings == 0);
	
	   if( nvars == -1 || vars == NULL )
	   {
	      int nbinvars;
	      int nintvars;
	      SCIP_CALL( SCIPgetVarsData(scip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	      nbinintvars = nbinvars + nintvars;
	      nvars = nbinintvars;
	   }
	   firstsol = sols[0];
	   assert(nvars > 0);
	
	   /* loop over integer and binary variables and check if their solution values match in all solutions */
	   for( v = 0; v < nvars; ++v )
	   {
	      SCIP_Real solval;
	      SCIP_VAR* var;
	      int s;
	
	      var = vars[v];
	      assert((v < SCIPgetNBinVars(scip) && SCIPvarIsBinary(var)) || (v >= SCIPgetNBinVars(scip) && SCIPvarIsIntegral(var)));
	      solval = SCIPgetSolVal(scip, firstsol, var);
	
	      /* determine if solution values match in all given solutions */
	      for( s = 1; s < nsols; ++s )
	      {
	         SCIP_Real solval2 = SCIPgetSolVal(scip, sols[s], var);
	         if( ! SCIPisEQ(scip, solval, solval2) )
	            break;
	      }
	
	      /* if we did not break early, all solutions agree on the solution value of this variable */
	      if( s == nsols )
	      {
	         tryAdd2variableBuffer(scip, var, solval, varbuf, valbuf, nfixings, TRUE);
	      }
	   }
	
	   return SCIP_OKAY;
	}
	
	/** callback to collect variable fixings of RINS */
	static
	DECL_VARFIXINGS(varFixingsRins)
	{
	   /*lint --e{715}*/
	   int nbinvars;
	   int nintvars;
	   SCIP_VAR** vars;
	   SCIP_SOL* incumbent;
	   SCIP_SOL* sols[2];
	   assert(scip != NULL);
	   assert(varbuf != NULL);
	   assert(nfixings != NULL);
	   assert(valbuf != NULL);
	
	   *result = SCIP_DELAYED;
	
	   if( ! SCIPhasCurrentNodeLP(scip) )
	      return SCIP_OKAY;
	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
	      return SCIP_OKAY;
	
	   *result = SCIP_DIDNOTRUN;
	
	   incumbent = SCIPgetBestSol(scip);
	   if( incumbent == NULL )
	      return SCIP_OKAY;
	
	   if( SCIPsolGetOrigin(incumbent) == SCIP_SOLORIGIN_ORIGINAL )
	      return SCIP_OKAY;
	
	   /* get variable information */
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	
	   /* return if no binary or integer variables are present */
	   if( nbinvars + nintvars == 0 )
	      return SCIP_OKAY;
	
	   sols[0] = NULL;
	   sols[1] = incumbent;
	
	   SCIP_CALL( fixMatchingSolutionValues(scip, sols, 2, vars, nbinvars + nintvars, varbuf, valbuf, nfixings) );
	
	   *result = SCIP_SUCCESS;
	
	   return SCIP_OKAY;
	}
	
	/** initialization callback for crossover when a new problem is read */
	static
	DECL_NHINIT(nhInitCrossover)
	{  /*lint --e{715}*/
	   DATA_CROSSOVER* data;
	
	   data = neighborhood->data.crossover;
	   assert(data != NULL);
	
	   if( data->rng != NULL )
	      SCIPfreeRandom(scip, &data->rng);
	
	   data->selsol = NULL;
	
	   SCIP_CALL( SCIPcreateRandom(scip, &data->rng, CROSSOVERSEED + (unsigned int)SCIPgetNVars(scip), TRUE) );
	
	   return SCIP_OKAY;
	}
	
	/** deinitialization callback for crossover when exiting a problem */
	static
	DECL_NHEXIT(nhExitCrossover)
	{  /*lint --e{715}*/
	   DATA_CROSSOVER* data;
	   data = neighborhood->data.crossover;
	
	   assert(neighborhood != NULL);
	   assert(data->rng != NULL);
	
	   SCIPfreeRandom(scip, &data->rng);
	
	   return SCIP_OKAY;
	}
	
	/** deinitialization callback for crossover before SCIP is freed */
	static
	DECL_NHFREE(nhFreeCrossover)
	{  /*lint --e{715}*/
	   assert(neighborhood->data.crossover != NULL);
	   SCIPfreeBlockMemory(scip, &neighborhood->data.crossover);
	
	   return SCIP_OKAY;
	}
	
	/** callback to collect variable fixings of crossover */
	static
	DECL_VARFIXINGS(varFixingsCrossover)
	{  /*lint --e{715}*/
	   DATA_CROSSOVER* data;
	   SCIP_RANDNUMGEN* rng;
	   SCIP_SOL** sols;
	   SCIP_SOL** scipsols;
	   int nsols;
	   int lastdraw;
	   assert(scip != NULL);
	   assert(varbuf != NULL);
	   assert(nfixings != NULL);
	   assert(valbuf != NULL);
	
	   data = neighborhood->data.crossover;
	
	   assert(data != NULL);
	   nsols = data->nsols;
	   data->selsol = NULL;
	
	   *result = SCIP_DIDNOTRUN;
	
	   /* return if the pool has not enough solutions */
	   if( nsols > SCIPgetNSols(scip) )
	      return SCIP_OKAY;
	
	   /* return if no binary or integer variables are present */
	   if( SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) == 0 )
	      return SCIP_OKAY;
	
	   rng = data->rng;
	   lastdraw = SCIPgetNSols(scip);
	   SCIP_CALL( SCIPallocBufferArray(scip, &sols, nsols) );
	   scipsols = SCIPgetSols(scip);
	
	   /* draw as many solutions from the pool as required by crossover, biased towards
	    * better solutions; therefore, the sorting of the solutions by objective is implicitly used
	    */
	   while( nsols > 0 )
	   {
	      /* no need for randomization anymore, exactly nsols many solutions remain for the selection */
	      if( lastdraw == nsols )
	      {
	         int s;
	
	         /* fill the remaining slots 0,...,nsols - 1 by the solutions at the same places */
	         for( s = 0; s < nsols; ++s )
	            sols[s] = scipsols[s];
	
	         nsols = 0;
	      }
	      else
	      {
	         int nextdraw;
	
	         assert(nsols < lastdraw);
	
	         /* draw from the lastdraw - nsols many solutions nsols - 1, ... lastdraw - 1 such that nsols many solution */
	         nextdraw = SCIPrandomGetInt(rng, nsols - 1, lastdraw - 1);
	         assert(nextdraw >= 0);
	
	         sols[nsols - 1] = scipsols[nextdraw];
	         nsols--;
	         lastdraw = nextdraw;
	      }
	   }
	
	   SCIP_CALL( fixMatchingSolutionValues(scip, sols, data->nsols, NULL, -1, varbuf, valbuf, nfixings) );
	
	   /* store best selected solution as reference solution */
	   data->selsol = sols[0];
	   assert(data->selsol != NULL);
	
	   *result = SCIP_SUCCESS;
	
	   SCIPfreeBufferArray(scip, &sols);
	
	   return SCIP_OKAY;
	}
	
	/** callback for crossover reference solution */
	static
	DECL_NHREFSOL(nhRefsolCrossover)
	{ /*lint --e{715}*/
	   DATA_CROSSOVER* data;
	
	   data = neighborhood->data.crossover;
	
	   if( data->selsol != NULL )
	   {
	      *solptr = data->selsol;
	      *result = SCIP_SUCCESS;
	   }
	   else
	   {
	      *result = SCIP_DIDNOTFIND;
	   }
	
	   return SCIP_OKAY;
	}
	
	/** initialization callback for mutation when a new problem is read */
	static
	DECL_NHINIT(nhInitMutation)
	{  /*lint --e{715}*/
	   DATA_MUTATION* data;
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	
	   SCIP_CALL( SCIPallocBlockMemory(scip, &neighborhood->data.mutation) );
	
	   data = neighborhood->data.mutation;
	   assert(data != NULL);
	
	   SCIP_CALL( SCIPcreateRandom(scip, &data->rng, MUTATIONSEED + (unsigned int)SCIPgetNVars(scip), TRUE) );
	
	   return SCIP_OKAY;
	}
	
	/** deinitialization callback for mutation when exiting a problem */
	static
	DECL_NHEXIT(nhExitMutation)
	{  /*lint --e{715}*/
	   DATA_MUTATION* data;
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	   data = neighborhood->data.mutation;
	   assert(data != NULL);
	
	   SCIPfreeRandom(scip, &data->rng);
	
	   SCIPfreeBlockMemory(scip, &neighborhood->data.mutation);
	
	   return SCIP_OKAY;
	}
	
	/** callback to collect variable fixings of mutation */
	static
	DECL_VARFIXINGS(varFixingsMutation)
	{  /*lint --e{715}*/
	   SCIP_RANDNUMGEN* rng;
	
	   SCIP_VAR** vars;
	   SCIP_VAR** varscpy;
	   int i;
	   int nvars;
	   int nbinvars;
	   int nintvars;
	   int nbinintvars;
	   int ntargetfixings;
	   SCIP_SOL* incumbentsol;
	   SCIP_Real targetfixingrate;
	
	   assert(scip != NULL);
	   assert(neighborhood != NULL);
	   assert(neighborhood->data.mutation != NULL);
	   assert(neighborhood->data.mutation->rng != NULL);
	   rng = neighborhood->data.mutation->rng;
	
	   *result = SCIP_DIDNOTRUN;
	
	   /* get the problem variables */
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, &nvars, &nbinvars, &nintvars, NULL, NULL) );
	
	   nbinintvars = nbinvars + nintvars;
	   if( nbinintvars == 0 )
	      return SCIP_OKAY;
	
	   incumbentsol = SCIPgetBestSol(scip);
	   if( incumbentsol == NULL )
	      return SCIP_OKAY;
	
	   targetfixingrate = neighborhood->fixingrate.targetfixingrate;
	   ntargetfixings = (int)(targetfixingrate * nbinintvars) + 1;
	
	   /* don't continue if number of discrete variables is too small to reach target fixing rate */
	   if( nbinintvars <= ntargetfixings )
	      return SCIP_OKAY;
	
	   *result = SCIP_DIDNOTFIND;
	
	   /* copy variables into a buffer array */
	   SCIP_CALL( SCIPduplicateBufferArray(scip, &varscpy, vars, nbinintvars) );
	
	   /* partially perturb the array until the number of target fixings is reached */
	   for( i = 0; *nfixings < ntargetfixings && i < nbinintvars; ++i )
	   {
	      int randint = SCIPrandomGetInt(rng, i, nbinintvars - 1);
	      assert(randint < nbinintvars);
	
	      if( randint > i )
	      {
	         SCIPswapPointers((void**)&varscpy[i], (void**)&varscpy[randint]);
	      }
	      /* copy the selected variables and their solution values into the buffer */
	      tryAdd2variableBuffer(scip, varscpy[i], SCIPgetSolVal(scip, incumbentsol, varscpy[i]), varbuf, valbuf, nfixings, TRUE);
	   }
	
	   assert(i == nbinintvars || *nfixings == ntargetfixings);
	
	   /* Not reaching the number of target fixings means that there is a significant fraction (at least 1 - targetfixingrate)
	    * of variables for which the incumbent solution value does not lie within the global bounds anymore. This is a nonsuccess
	    * for the neighborhood (additional fixings are not possible), which is okay because the incumbent solution is
	    * significantly outdated
	    */
	   if( *nfixings == ntargetfixings )
	      *result = SCIP_SUCCESS;
	
	   /* free the buffer array */
	   SCIPfreeBufferArray(scip, &varscpy);
	
	   return SCIP_OKAY;
	}
	
	/** add local branching constraint */
	static
	SCIP_RETCODE addLocalBranchingConstraint(
	   SCIP*                 sourcescip,         /**< source SCIP data structure */
	   SCIP*                 targetscip,         /**< target SCIP data structure */
	   SCIP_VAR**            subvars,            /**< array of sub SCIP variables in same order as source SCIP variables */
	   int                   distance,           /**< right hand side of the local branching constraint */
	   SCIP_Bool*            success,            /**< pointer to store of a local branching constraint has been successfully added */
	   int*                  naddedconss         /**< pointer to increase the number of added constraints */
	   )
	{
	   int nbinvars;
	   int i;
	   SCIP_SOL* referencesol;
	   SCIP_CONS* localbranchcons;
	   SCIP_VAR** vars;
	   SCIP_Real* consvals;
	   SCIP_Real rhs;
	
	   assert(sourcescip != NULL);
	   assert(*success == FALSE);
	
	   nbinvars = SCIPgetNBinVars(sourcescip);
	   vars = SCIPgetVars(sourcescip);
	
	   if( nbinvars <= 3 )
	      return SCIP_OKAY;
	
	   referencesol = SCIPgetBestSol(sourcescip);
	   if( referencesol == NULL )
	      return SCIP_OKAY;
	
	   rhs = (SCIP_Real)distance;
	   rhs = MAX(rhs, 2.0);
	
	   SCIP_CALL( SCIPallocBufferArray(sourcescip, &consvals, nbinvars) );
	
	   /* loop over binary variables and fill the local branching constraint */
	   for( i = 0; i < nbinvars; ++i )
	   {
	      /* skip variables that are not present in sub-SCIP */
	      if( subvars[i] == NULL )
	         continue;
	
	      if( SCIPisEQ(sourcescip, SCIPgetSolVal(sourcescip, referencesol, vars[i]), 0.0) )
	         consvals[i] = 1.0;
	      else
	      {
	         consvals[i] = -1.0;
	         rhs -= 1.0;
	      }
	   }
	
	   /* create the local branching constraint in the target scip */
	   SCIP_CALL( SCIPcreateConsBasicLinear(targetscip, &localbranchcons, "localbranch", nbinvars, subvars, consvals, -SCIPinfinity(sourcescip), rhs) );
	   SCIP_CALL( SCIPaddCons(targetscip, localbranchcons) );
	   SCIP_CALL( SCIPreleaseCons(targetscip, &localbranchcons) );
	
	   *naddedconss = 1;
	   *success = TRUE;
	
	   SCIPfreeBufferArray(sourcescip, &consvals);
	
	   return SCIP_OKAY;
	}
	
	/** callback for local branching subproblem changes */
	static
	DECL_CHANGESUBSCIP(changeSubscipLocalbranching)
	{  /*lint --e{715}*/
	
	   SCIP_CALL( addLocalBranchingConstraint(sourcescip, targetscip, subvars, (int)(0.2 * SCIPgetNBinVars(sourcescip)), success, naddedconss) );
	
	   return SCIP_OKAY;
	}
	
	/** callback for proximity subproblem changes */
	static
	DECL_CHANGESUBSCIP(changeSubscipProximity)
	{  /*lint --e{715}*/
	   SCIP_SOL* referencesol;
	   SCIP_VAR** vars;
	   int nbinvars;
	   int nintvars;
	   int nvars;
	   int i;
	
	   SCIP_CALL( SCIPgetVarsData(sourcescip, &vars, &nvars, &nbinvars, &nintvars, NULL, NULL) );
	
	   if( nbinvars == 0 )
	      return SCIP_OKAY;
	
	   referencesol = SCIPgetBestSol(sourcescip);
	   if( referencesol == NULL )
	      return SCIP_OKAY;
	
	   /* loop over binary variables, set objective coefficients based on reference solution in a local branching fashion */
	   for( i = 0; i < nbinvars; ++i )
	   {
	      SCIP_Real newobj;
	
	      /* skip variables not present in sub-SCIP */
	      if( subvars[i] == NULL )
	         continue;
	
	      if( SCIPgetSolVal(sourcescip, referencesol, vars[i]) < 0.5 )
	         newobj = -1.0;
	      else
	         newobj = 1.0;
	      SCIP_CALL( SCIPchgVarObj(targetscip, subvars[i], newobj) );
	   }
	
	   /* loop over the remaining variables and change their objective coefficients to 0 */
	   for( ; i < nvars; ++i )
	   {
	      /* skip variables not present in sub-SCIP */
	      if( subvars[i] == NULL )
	         continue;
	
	      SCIP_CALL( SCIPchgVarObj(targetscip, subvars[i], 0.0) );
	   }
	
	   *nchgobjs = nvars;
	   *success = TRUE;
	
	   return SCIP_OKAY;
	}
	
	/** callback for zeroobjective subproblem changes */
	static
	DECL_CHANGESUBSCIP(changeSubscipZeroobjective)
	{  /*lint --e{715}*/
	   SCIP_CONSHDLR* conshdlrnl;
	   SCIP_VAR** vars;
	   int nvars;
	   int i;
	
	   assert(*success == FALSE);
	
	   SCIP_CALL( SCIPgetVarsData(sourcescip, &vars, &nvars, NULL, NULL, NULL, NULL) );
	
	   /* do not run if no objective variables are present */
	   if( SCIPgetNObjVars(sourcescip) == 0 )
	      return SCIP_OKAY;
	
	   /* zeroobj may trigger fixing objvar in nonlinear constraint to infinity,
	    * which expr_var.c:simplify cannot handle at the moment; also #3273
	    */
	   conshdlrnl = SCIPfindConshdlr(sourcescip, "nonlinear");
	   if( conshdlrnl != NULL && SCIPconshdlrGetNActiveConss(conshdlrnl) > 0 )
	      return SCIP_OKAY;
	
	   /* loop over the variables and change their objective coefficients to 0 */
	   for( i = 0; i < nvars; ++i )
	   {
	      /* skip variables not present in sub-SCIP */
	      if( subvars[i] == NULL )
	         continue;
	
	      SCIP_CALL( SCIPchgVarObj(targetscip, subvars[i], 0.0) );
	   }
	
	   *nchgobjs = nvars;
	   *success = TRUE;
	
	   return SCIP_OKAY;
	}
	
	/** compute tightened bounds for integer variables depending on how much the LP and the incumbent solution values differ */
	static
	void computeIntegerVariableBoundsDins(
	   SCIP*                 scip,               /**< SCIP data structure of the original problem */
	   SCIP_VAR*             var,                /**< the variable for which bounds should be computed */
	   SCIP_Real*            lbptr,              /**< pointer to store the lower bound in the DINS sub-SCIP */
	   SCIP_Real*            ubptr               /**< pointer to store the upper bound in the DINS sub-SCIP */
	   )
	{
	   SCIP_Real mipsol;
	   SCIP_Real lpsol;
	
	   SCIP_Real lbglobal;
	   SCIP_Real ubglobal;
	   SCIP_SOL* bestsol;
	
	   /* get the bounds for each variable */
	   lbglobal = SCIPvarGetLbGlobal(var);
	   ubglobal = SCIPvarGetUbGlobal(var);
	
	   assert(SCIPvarGetType(var) == SCIP_VARTYPE_INTEGER);
	   /* get the current LP solution for each variable */
	   lpsol = SCIPvarGetLPSol(var);
	
	   /* get the current MIP solution for each variable */
	   bestsol = SCIPgetBestSol(scip);
	   mipsol = SCIPgetSolVal(scip, bestsol, var);
	
	   /* if the solution values differ by 0.5 or more, the variable is rebounded, otherwise it is just copied */
	   if( REALABS(lpsol - mipsol) >= 0.5 )
	   {
	      SCIP_Real range;
	
	      *lbptr = lbglobal;
	      *ubptr = ubglobal;
	
	      /* create an equally sized range around lpsol for general integers: bounds are lpsol +- (mipsol-lpsol) */
	      range = 2 * lpsol - mipsol;
	
	      if( mipsol >= lpsol )
	      {
	         range = SCIPfeasCeil(scip, range);
	         *lbptr = MAX(*lbptr, range);
	
	         /* when the bound new upper bound is equal to the current MIP solution, we set both bounds to the integral bound (without eps) */
	         if( SCIPisFeasEQ(scip, mipsol, *lbptr) )
	            *ubptr = *lbptr;
	         else
	            *ubptr = mipsol;
	      }
	      else
	      {
	         range = SCIPfeasFloor(scip, range);
	         *ubptr = MIN(*ubptr, range);
	
	         /* when the bound new upper bound is equal to the current MIP solution, we set both bounds to the integral bound (without eps) */
	         if( SCIPisFeasEQ(scip, mipsol, *ubptr) )
	            *lbptr = *ubptr;
	         else
	            *lbptr = mipsol;
	      }
	
	      /* the global domain of variables might have been reduced since incumbent was found: adjust lb and ub accordingly */
	      *lbptr = MAX(*lbptr, lbglobal);
	      *ubptr = MIN(*ubptr, ubglobal);
	   }
	   else
	   {
	      /* the global domain of variables might have been reduced since incumbent was found: adjust it accordingly */
	      *lbptr = MAX(mipsol, lbglobal);
	      *ubptr = MIN(mipsol, ubglobal);
	   }
	}
	
	/** callback to collect variable fixings of DINS */
	static
	DECL_VARFIXINGS(varFixingsDins)
	{
	   DATA_DINS* data;
	   SCIP_SOL* rootlpsol;
	   SCIP_SOL** sols;
	   int nsols;
	   int nmipsols;
	   int nbinvars;
	   int nintvars;
	   SCIP_VAR** vars;
	   int v;
	
	   data = neighborhood->data.dins;
	   assert(data != NULL);
	   nmipsols = SCIPgetNSols(scip);
	   nmipsols = MIN(nmipsols, data->npoolsols);
	
	   *result = SCIP_DELAYED;
	
	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
	      return SCIP_OKAY;
	
	   *result = SCIP_DIDNOTRUN;
	
	   if( nmipsols == 0 )
	      return SCIP_OKAY;
	
	   SCIP_CALL( SCIPgetVarsData(scip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	
	   if( nbinvars + nintvars == 0 )
	      return SCIP_OKAY;
	
	   SCIP_CALL( SCIPcreateSol(scip, &rootlpsol, NULL) );
	
	   /* save root solution LP values in solution */
	   for( v = 0; v < nbinvars + nintvars; ++v )
	   {
	      SCIP_CALL( SCIPsetSolVal(scip, rootlpsol, vars[v], SCIPvarGetRootSol(vars[v])) );
	   }
	
	   /* add the node and the root LP solution */
	   nsols = nmipsols + 2;
	
	   SCIP_CALL( SCIPallocBufferArray(scip, &sols, nsols) );
	   sols[0] = NULL; /* node LP solution */
	   sols[1] = rootlpsol;
	
	   /* copy the remaining MIP solutions after the LP solutions */
	   BMScopyMemoryArray(&sols[2], SCIPgetSols(scip), nmipsols); /*lint !e866*/
	
	   /* 1. Binary variables are fixed if their values agree in all the solutions */
	   if( nbinvars > 0 )
	   {
	      SCIP_CALL( fixMatchingSolutionValues(scip, sols, nsols, vars, nbinvars, varbuf, valbuf, nfixings) );
	   }
	
	   /* 2. Integer variables are fixed if they have a very low distance between the incumbent and the root LP solution */
	   for( v = nbinvars; v < nintvars; ++v )
	   {
	      SCIP_Real lb;
	      SCIP_Real ub;
	      computeIntegerVariableBoundsDins(scip, vars[v], &lb, &ub);
	
	      if( ub - lb < 0.5 )
	      {
	         assert(SCIPisFeasIntegral(scip, lb));
	         tryAdd2variableBuffer(scip, vars[v], lb, varbuf, valbuf, nfixings, TRUE);
	      }
	   }
	
	   *result = SCIP_SUCCESS;
	
	   SCIPfreeBufferArray(scip, &sols);
	
	   SCIP_CALL( SCIPfreeSol(scip, &rootlpsol) );
	
	   return SCIP_OKAY;
	}
	
	/** callback for DINS subproblem changes */
	static
	DECL_CHANGESUBSCIP(changeSubscipDins)
	{  /*lint --e{715}*/
	   SCIP_VAR** vars;
	   int nintvars;
	   int nbinvars;
	   int v;
	
	   SCIP_CALL( SCIPgetVarsData(sourcescip, &vars, NULL, &nbinvars, &nintvars, NULL, NULL) );
	
	   /* 1. loop over integer variables and tighten the bounds */
	   for( v = nbinvars; v < nintvars; ++v )
	   {
	      SCIP_Real lb;
	      SCIP_Real ub;
	
	      /* skip variables not present in sub-SCIP */
	      if( subvars[v] == NULL )
	         continue;
	
	      computeIntegerVariableBoundsDins(sourcescip, vars[v], &lb, &ub);
	
	      SCIP_CALL( SCIPchgVarLbGlobal(targetscip, subvars[v], lb) );
	      SCIP_CALL( SCIPchgVarUbGlobal(targetscip, subvars[v], ub) );
	      ++(*ndomchgs);
	   }
	
	   /* 2. add local branching constraint for binary variables */
	   SCIP_CALL( addLocalBranchingConstraint(sourcescip, targetscip, subvars, (int)(0.1 * SCIPgetNBinVars(sourcescip)), success, naddedconss) );
	
	   *success = TRUE;
	
	   return SCIP_OKAY;
	}
	
	/** deinitialization callback for DINS before SCIP is freed */
	static
	DECL_NHFREE(nhFreeDins)
	{
	   assert(neighborhood->data.dins != NULL);
	
	   SCIPfreeBlockMemory(scip, &neighborhood->data.dins);
	
	   return SCIP_OKAY;
	}
	
	/** deinitialization callback for trustregion before SCIP is freed */
	static
	DECL_NHFREE(nhFreeTrustregion)
	{
	   assert(neighborhood->data.trustregion != NULL);
	
	   SCIPfreeBlockMemory(scip, &neighborhood->data.trustregion);
	
	   return SCIP_OKAY;
	}
	
	/** add trust region neighborhood constraint and auxiliary objective variable */
	static
	DECL_CHANGESUBSCIP(changeSubscipTrustregion)
	{  /*lint --e{715}*/
	   DATA_TRUSTREGION* data;
	
	   data = neighborhood->data.trustregion;
	
	   /* adding the neighborhood constraint for the trust region heuristic */
	   SCIP_CALL( SCIPaddTrustregionNeighborhoodConstraint(sourcescip, targetscip, subvars, data->violpenalty) );
	
	   /* incrementing the change in objective since an additional variable is added to the objective to penalize the
	    * violation of the trust region.
	    */
	   ++(*nchgobjs);
	
	   return SCIP_OKAY;
	}
	
	/** callback that returns the incumbent solution as a reference point */
	static
	DECL_NHREFSOL(nhRefsolIncumbent)
	{  /*lint --e{715}*/
	   assert(scip != NULL);
	
	   if( SCIPgetBestSol(scip) != NULL )
	   {
	      *result = SCIP_SUCCESS;
	      *solptr = SCIPgetBestSol(scip);
	   }
	   else
	   {
	      *result = SCIP_DIDNOTFIND;
	   }
	
	   return SCIP_OKAY;
	}
	
	
	/** callback function that deactivates a neighborhood on problems with no discrete variables */
	static
	DECL_NHDEACTIVATE(nhDeactivateDiscreteVars)
	{ /*lint --e{715}*/
	   assert(scip != NULL);
	   assert(deactivate != NULL);
	
	   /* deactivate if no discrete variables are present */
	   *deactivate = (SCIPgetNBinVars(scip) + SCIPgetNIntVars(scip) == 0);
	
	   return SCIP_OKAY;
	}
	
	/** callback function that deactivates a neighborhood on problems with no binary variables */
	static
	DECL_NHDEACTIVATE(nhDeactivateBinVars)
	{ /*lint --e{715}*/
	   assert(scip != NULL);
	   assert(deactivate != NULL);
	
	   /* deactivate if no discrete variables are present */
	   *deactivate = (SCIPgetNBinVars(scip) == 0);
	
	   return SCIP_OKAY;
	}
	
	/** callback function that deactivates a neighborhood on problems with no objective variables */
	static
	DECL_NHDEACTIVATE(nhDeactivateObjVars)
	{ /*lint --e{715}*/
	   assert(scip != NULL);
	   assert(deactivate != NULL);
	
	   /* deactivate if no discrete variables are present */
	   *deactivate = (SCIPgetNObjVars(scip) == 0);
	
	   return SCIP_OKAY;
	}
	
	
	/** include all neighborhoods */
	static
	SCIP_RETCODE includeNeighborhoods(
	   SCIP*                 scip,               /**< SCIP data structure */
	   SCIP_HEURDATA*        heurdata            /**< heuristic data of the ALNS heuristic */
	   )
	{
	   NH* rens;
	   NH* rins;
	   NH* mutation;
	   NH* localbranching;
	   NH* crossover;
	   NH* proximity;
	   NH* zeroobjective;
	   NH* dins;
	   NH* trustregion;
	
	   heurdata->nneighborhoods = 0;
	
	   /* include the RENS neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &rens, "rens",
	         DEFAULT_MINFIXINGRATE_RENS, DEFAULT_MAXFIXINGRATE_RENS, DEFAULT_ACTIVE_RENS, DEFAULT_PRIORITY_RENS,
	         varFixingsRens, changeSubscipRens, NULL, NULL, NULL, NULL, nhDeactivateDiscreteVars) );
	
	   /* include the RINS neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &rins, "rins",
	         DEFAULT_MINFIXINGRATE_RINS, DEFAULT_MAXFIXINGRATE_RINS, DEFAULT_ACTIVE_RINS, DEFAULT_PRIORITY_RINS,
	         varFixingsRins, NULL, NULL, NULL, NULL, nhRefsolIncumbent, nhDeactivateDiscreteVars) );
	
	   /* include the mutation neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &mutation, "mutation",
	         DEFAULT_MINFIXINGRATE_MUTATION, DEFAULT_MAXFIXINGRATE_MUTATION, DEFAULT_ACTIVE_MUTATION, DEFAULT_PRIORITY_MUTATION,
	         varFixingsMutation, NULL, nhInitMutation, nhExitMutation, NULL, nhRefsolIncumbent, nhDeactivateDiscreteVars) );
	
	   /* include the local branching neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &localbranching, "localbranching",
	         DEFAULT_MINFIXINGRATE_LOCALBRANCHING, DEFAULT_MAXFIXINGRATE_LOCALBRANCHING, DEFAULT_ACTIVE_LOCALBRANCHING, DEFAULT_PRIORITY_LOCALBRANCHING,
	         NULL, changeSubscipLocalbranching, NULL, NULL, NULL, nhRefsolIncumbent, nhDeactivateBinVars) );
	
	   /* include the crossover neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &crossover, "crossover",
	         DEFAULT_MINFIXINGRATE_CROSSOVER, DEFAULT_MAXFIXINGRATE_CROSSOVER, DEFAULT_ACTIVE_CROSSOVER, DEFAULT_PRIORITY_CROSSOVER,
	         varFixingsCrossover, NULL,
	         nhInitCrossover, nhExitCrossover, nhFreeCrossover, nhRefsolCrossover, nhDeactivateDiscreteVars) );
	
	   /* allocate data for crossover to include the parameter */
	   SCIP_CALL( SCIPallocBlockMemory(scip, &crossover->data.crossover) );
	   crossover->data.crossover->rng = NULL;
	
	   /* add crossover neighborhood parameters */
	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/alns/crossover/nsols", "the number of solutions that crossover should combine",
	         &crossover->data.crossover->nsols, TRUE, DEFAULT_NSOLS_CROSSOVER, 2, 10, NULL, NULL) );
	
	   /* include the Proximity neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &proximity, "proximity",
	         DEFAULT_MINFIXINGRATE_PROXIMITY, DEFAULT_MAXFIXINGRATE_PROXIMITY, DEFAULT_ACTIVE_PROXIMITY, DEFAULT_PRIORITY_PROXIMITY,
	         NULL, changeSubscipProximity, NULL, NULL, NULL, nhRefsolIncumbent, nhDeactivateBinVars) );
	
	   /* include the Zeroobjective neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &zeroobjective, "zeroobjective",
	         DEFAULT_MINFIXINGRATE_ZEROOBJECTIVE, DEFAULT_MAXFIXINGRATE_ZEROOBJECTIVE, DEFAULT_ACTIVE_ZEROOBJECTIVE, DEFAULT_PRIORITY_ZEROOBJECTIVE,
	         NULL, changeSubscipZeroobjective, NULL, NULL, NULL, nhRefsolIncumbent, nhDeactivateObjVars) );
	
	   /* include the DINS neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &dins, "dins",
	         DEFAULT_MINFIXINGRATE_DINS, DEFAULT_MAXFIXINGRATE_DINS, DEFAULT_ACTIVE_DINS, DEFAULT_PRIORITY_DINS,
	         varFixingsDins, changeSubscipDins, NULL, NULL, nhFreeDins, nhRefsolIncumbent, nhDeactivateBinVars) );
	
	   /* allocate data for DINS to include the parameter */
	   SCIP_CALL( SCIPallocBlockMemory(scip, &dins->data.dins) );
	
	   /* add DINS neighborhood parameters */
	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/alns/dins/npoolsols",
	         "number of pool solutions where binary solution values must agree",
	         &dins->data.dins->npoolsols, TRUE, DEFAULT_NPOOLSOLS_DINS, 1, 100, NULL, NULL) );
	
	   /* include the trustregion neighborhood */
	   SCIP_CALL( alnsIncludeNeighborhood(scip, heurdata, &trustregion, "trustregion",
	         DEFAULT_MINFIXINGRATE_TRUSTREGION, DEFAULT_MAXFIXINGRATE_TRUSTREGION, DEFAULT_ACTIVE_TRUSTREGION, DEFAULT_PRIORITY_TRUSTREGION,
	         NULL, changeSubscipTrustregion, NULL, NULL, nhFreeTrustregion, nhRefsolIncumbent, nhDeactivateBinVars) );
	
	   /* allocate data for trustregion to include the parameter */
	   SCIP_CALL( SCIPallocBlockMemory(scip, &trustregion->data.trustregion) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/trustregion/violpenalty",
	         "the penalty for each change in the binary variables from the candidate solution",
	         &trustregion->data.trustregion->violpenalty, FALSE, DEFAULT_VIOLPENALTY_TRUSTREGION, 0.0, SCIP_REAL_MAX, NULL, NULL) );
	
	   return SCIP_OKAY;
	}
	
	/** initialization method of primal heuristic (called after problem was transformed) */
	static
	SCIP_DECL_HEURINIT(heurInitAlns)
	{  /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	   int i;
	
	   assert(scip != NULL);
	   assert(heur != NULL);
	
	   heurdata = SCIPheurGetData(heur);
	   assert(heurdata != NULL);
	
	   /* reactivate all neighborhoods if a new problem is read in */
	   heurdata->nactiveneighborhoods = heurdata->nneighborhoods;
	
	   /* initialize neighborhoods for new problem */
	   for( i = 0; i < heurdata->nneighborhoods; ++i )
	   {
	      NH* neighborhood = heurdata->neighborhoods[i];
	
	      SCIP_CALL( neighborhoodInit(scip, neighborhood) );
	
	      SCIP_CALL( resetFixingRate(scip, &neighborhood->fixingrate) );
	
	      SCIP_CALL( neighborhoodStatsReset(scip, &neighborhood->stats) );
	   }
	
	   /* open reward file for reading */
	   if( strncasecmp(heurdata->rewardfilename, DEFAULT_REWARDFILENAME, strlen(DEFAULT_REWARDFILENAME)) != 0 )
	   {
	      heurdata->rewardfile = fopen(heurdata->rewardfilename, "w");
	
	      if( heurdata->rewardfile == NULL )
	      {
	         SCIPerrorMessage("Error: Could not open reward file <%s>\n", heurdata->rewardfilename);
	         return SCIP_FILECREATEERROR;
	      }
	
	      SCIPdebugMsg(scip, "Writing reward information to <%s>\n", heurdata->rewardfilename);
	   }
	   else
	      heurdata->rewardfile = NULL;
	
	   return SCIP_OKAY;
	}
	
	
	/** solving process initialization method of primal heuristic (called when branch and bound process is about to begin) */
	static
	SCIP_DECL_HEURINITSOL(heurInitsolAlns)
	{  /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	   int i;
	   SCIP_Real* priorities;
	   unsigned int initseed;
	
	   assert(scip != NULL);
	   assert(heur != NULL);
	
	   heurdata = SCIPheurGetData(heur);
	   assert(heurdata != NULL);
	   heurdata->nactiveneighborhoods = heurdata->nneighborhoods;
	
	   SCIP_CALL( SCIPallocBufferArray(scip, &priorities, heurdata->nactiveneighborhoods) );
	
	   /* init neighborhoods for new problem by resetting their statistics and fixing rate */
	   for( i = heurdata->nneighborhoods - 1; i >= 0; --i )
	   {
	      NH* neighborhood = heurdata->neighborhoods[i];
	      SCIP_Bool deactivate;
	
	      SCIP_CALL( neighborhood->nhdeactivate(scip, &deactivate) );
	
	      /* disable inactive neighborhoods */
	      if( deactivate || ! neighborhood->active )
	      {
	         if( heurdata->nactiveneighborhoods - 1 > i )
	         {
	            assert(heurdata->neighborhoods[heurdata->nactiveneighborhoods - 1]->active);
	            SCIPswapPointers((void **)&heurdata->neighborhoods[i], (void **)&heurdata->neighborhoods[heurdata->nactiveneighborhoods - 1]);
	         }
	         heurdata->nactiveneighborhoods--;
	      }
	   }
	
	   /* collect neighborhood priorities */
	   for( i = 0; i < heurdata->nactiveneighborhoods; ++i )
	      priorities[i] = heurdata->neighborhoods[i]->priority;
	
	   initseed = (unsigned int)(heurdata->seed + SCIPgetNVars(scip));
	
	   /* active neighborhoods might change between init calls, reset functionality must take this into account */
	   if( heurdata->bandit != NULL && SCIPbanditGetNActions(heurdata->bandit) != heurdata->nactiveneighborhoods )
	   {
	      SCIP_CALL( SCIPfreeBandit(scip, &heurdata->bandit) );
	
	      heurdata->bandit = NULL;
	   }
	
	   if( heurdata->nactiveneighborhoods > 0 )
	   {  /* create or reset bandit algorithm */
	      if( heurdata->bandit == NULL )
	      {
	         SCIP_CALL( createBandit(scip, heurdata, priorities, initseed) );
	
	         resetMinimumImprovement(heurdata);
	         resetTargetNodeLimit(heurdata);
	      }
	      else if( heurdata->resetweights )
	      {
	         SCIP_CALL( SCIPresetBandit(scip, heurdata->bandit, priorities, initseed) );
	
	         resetMinimumImprovement(heurdata);
	         resetTargetNodeLimit(heurdata);
	      }
	   }
	
	   heurdata->usednodes = 0;
	   heurdata->ninitneighborhoods = heurdata->nactiveneighborhoods;
	
	   heurdata->lastcallsol = NULL;
	   heurdata->firstcallthissol = 0;
	
	   resetCurrentNeighborhood(heurdata);
	
	   SCIPfreeBufferArray(scip, &priorities);
	
	   return SCIP_OKAY;
	}
	
	
	/** deinitialization method of primal heuristic (called before transformed problem is freed) */
	static
	SCIP_DECL_HEUREXIT(heurExitAlns)
	{  /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	   int i;
	
	   assert(scip != NULL);
	   assert(heur != NULL);
	
	   heurdata = SCIPheurGetData(heur);
	   assert(heurdata != NULL);
	
	   /* free neighborhood specific data */
	   for( i = 0; i < heurdata->nneighborhoods; ++i )
	   {
	      NH* neighborhood = heurdata->neighborhoods[i];
	
	      SCIP_CALL( neighborhoodExit(scip, neighborhood) );
	   }
	
	   if( heurdata->rewardfile != NULL )
	   {
	      fclose(heurdata->rewardfile);
	      heurdata->rewardfile = NULL;
	   }
	
	   return SCIP_OKAY;
	}
	
	/** destructor of primal heuristic to free user data (called when SCIP is exiting) */
	static
	SCIP_DECL_HEURFREE(heurFreeAlns)
	{  /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	   int i;
	
	   assert(scip != NULL);
	   assert(heur != NULL);
	
	   heurdata = SCIPheurGetData(heur);
	   assert(heurdata != NULL);
	
	   /* bandits are only initialized if a problem has been read */
	   if( heurdata->bandit != NULL )
	   {
	      SCIP_CALL( SCIPfreeBandit(scip, &heurdata->bandit) );
	   }
	
	   /* free neighborhoods */
	   for( i = 0; i < heurdata->nneighborhoods; ++i )
	   {
	      SCIP_CALL( alnsFreeNeighborhood(scip, &(heurdata->neighborhoods[i])) );
	   }
	
	   SCIPfreeBlockMemoryArray(scip, &heurdata->neighborhoods, NNEIGHBORHOODS);
	
	   SCIPfreeBlockMemory(scip, &heurdata);
	
	   return SCIP_OKAY;
	}
	
	/** output method of statistics table to output file stream 'file' */
	static
	SCIP_DECL_TABLEOUTPUT(tableOutputNeighborhood)
	{ /*lint --e{715}*/
	   SCIP_HEURDATA* heurdata;
	
	   assert(SCIPfindHeur(scip, HEUR_NAME) != NULL);
	   heurdata = SCIPheurGetData(SCIPfindHeur(scip, HEUR_NAME));
	   assert(heurdata != NULL);
	
	   printNeighborhoodStatistics(scip, heurdata, file);
	
	   return SCIP_OKAY;
	}
	
	/*
	 * primal heuristic specific interface methods
	 */
	
	/** creates the alns primal heuristic and includes it in SCIP */
	SCIP_RETCODE SCIPincludeHeurAlns(
	   SCIP*                 scip                /**< SCIP data structure */
	   )
	{
	   SCIP_HEURDATA* heurdata;
	   SCIP_HEUR* heur;
	
	   /* create alns primal heuristic data */
	   heurdata = NULL;
	   heur = NULL;
	
	   SCIP_CALL( SCIPallocBlockMemory(scip, &heurdata) );
	   BMSclearMemory(heurdata);
	
	   /* TODO make this a user parameter? */
	   heurdata->lplimfac = LPLIMFAC;
	
	   SCIP_CALL( SCIPallocBlockMemoryArray(scip, &heurdata->neighborhoods, NNEIGHBORHOODS) );
	
	   /* include primal heuristic */
	   SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
	         HEUR_NAME, HEUR_DESC, HEUR_DISPCHAR, HEUR_PRIORITY, HEUR_FREQ, HEUR_FREQOFS,
	         HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecAlns, heurdata) );
	
	   assert(heur != NULL);
	
	   /* include all neighborhoods */
	   SCIP_CALL( includeNeighborhoods(scip, heurdata) );
	
	   /* set non fundamental callbacks via setter functions */
	   SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyAlns) );
	   SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeAlns) );
	   SCIP_CALL( SCIPsetHeurInit(scip, heur, heurInitAlns) );
	   SCIP_CALL( SCIPsetHeurInitsol(scip, heur, heurInitsolAlns) );
	   SCIP_CALL( SCIPsetHeurExit(scip, heur, heurExitAlns) );
	
	   /* add alns primal heuristic parameters */
	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/maxnodes",
	         "maximum number of nodes to regard in the subproblem",
	         &heurdata->maxnodes,  TRUE,DEFAULT_MAXNODES, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/nodesofs",
	         "offset added to the nodes budget",
	         &heurdata->nodesoffset, FALSE, DEFAULT_NODESOFFSET, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/minnodes",
	         "minimum number of nodes required to start a sub-SCIP",
	         &heurdata->minnodes, TRUE, DEFAULT_MINNODES, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddLongintParam(scip, "heuristics/" HEUR_NAME "/waitingnodes",
	         "number of nodes since last incumbent solution that the heuristic should wait",
	         &heurdata->waitingnodes, TRUE, DEFAULT_WAITINGNODES, 0LL, SCIP_LONGINT_MAX, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/nodesquot",
	         "fraction of nodes compared to the main SCIP for budget computation",
	         &heurdata->nodesquot, FALSE, DEFAULT_NODESQUOT, 0.0, 1.0, NULL, NULL) );
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/nodesquotmin",
	         "lower bound fraction of nodes compared to the main SCIP for budget computation",
	         &heurdata->nodesquotmin, FALSE, DEFAULT_NODESQUOTMIN, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/startminimprove",
	         "initial factor by which ALNS should at least improve the incumbent",
	         &heurdata->startminimprove, TRUE, DEFAULT_STARTMINIMPROVE, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/minimprovelow",
	         "lower threshold for the minimal improvement over the incumbent",
	         &heurdata->minimprovelow, TRUE, DEFAULT_MINIMPROVELOW, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/minimprovehigh",
	         "upper bound for the minimal improvement over the incumbent",
	         &heurdata->minimprovehigh, TRUE, DEFAULT_MINIMPROVEHIGH, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/nsolslim",
	         "limit on the number of improving solutions in a sub-SCIP call",
	         &heurdata->nsolslim, FALSE, DEFAULT_NSOLSLIM, -1, INT_MAX, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddCharParam(scip, "heuristics/" HEUR_NAME "/banditalgo",
	         "the bandit algorithm: (u)pper confidence bounds, (e)xp.3, epsilon (g)reedy",
	         &heurdata->banditalgo, TRUE, DEFAULT_BANDITALGO, "ueg", NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/gamma",
	         "weight between uniform (gamma ~ 1) and weight driven (gamma ~ 0) probability distribution for exp3",
	         &heurdata->exp3_gamma, TRUE, DEFAULT_GAMMA, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/beta",
	         "reward offset between 0 and 1 at every observation for Exp.3",
	         &heurdata->exp3_beta, TRUE, DEFAULT_BETA, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/alpha",
	         "parameter to increase the confidence width in UCB",
	         &heurdata->ucb_alpha, TRUE, DEFAULT_ALPHA, 0.0, 100.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/usedistances",
	         "distances from fixed variables be used for variable prioritization",
	         &heurdata->usedistances, TRUE, DEFAULT_USEDISTANCES, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/useredcost",
	         "should reduced cost scores be used for variable prioritization?",
	         &heurdata->useredcost, TRUE, DEFAULT_USEREDCOST, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/domorefixings",
	         "should the ALNS heuristic do more fixings by itself based on variable prioritization "
	         "until the target fixing rate is reached?",
	         &heurdata->domorefixings, TRUE, DEFAULT_DOMOREFIXINGS, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/adjustfixingrate",
	         "should the heuristic adjust the target fixing rate based on the success?",
	         &heurdata->adjustfixingrate, TRUE, DEFAULT_ADJUSTFIXINGRATE, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/usesubscipheurs",
	         "should the heuristic activate other sub-SCIP heuristics during its search?",
	         &heurdata->usesubscipheurs, TRUE, DEFAULT_USESUBSCIPHEURS, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/rewardcontrol",
	         "reward control to increase the weight of the simple solution indicator and decrease the weight of the closed gap reward",
	         &heurdata->rewardcontrol, TRUE, DEFAULT_REWARDCONTROL, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/targetnodefactor",
	         "factor by which target node number is eventually increased",
	         &heurdata->targetnodefactor, TRUE, DEFAULT_TARGETNODEFACTOR, 1.0, 1e+5, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/seed",
	         "initial random seed for bandit algorithms and random decisions by neighborhoods",
	         &heurdata->seed, FALSE, DEFAULT_SEED, 0, INT_MAX, NULL, NULL) );
	   SCIP_CALL( SCIPaddIntParam(scip, "heuristics/" HEUR_NAME "/maxcallssamesol",
	         "number of allowed executions of the heuristic on the same incumbent solution (-1: no limit, 0: number of active neighborhoods)",
	         &heurdata->maxcallssamesol, TRUE, DEFAULT_MAXCALLSSAMESOL, -1, 100, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/adjustminimprove",
	         "should the factor by which the minimum improvement is bound be dynamically updated?",
	         &heurdata->adjustminimprove, TRUE, DEFAULT_ADJUSTMINIMPROVE, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/adjusttargetnodes",
	         "should the target nodes be dynamically adjusted?",
	         &heurdata->adjusttargetnodes, TRUE, DEFAULT_ADJUSTTARGETNODES, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/eps",
	         "increase exploration in epsilon-greedy bandit algorithm",
	         &heurdata->epsgreedy_eps, TRUE, DEFAULT_EPS, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/rewardbaseline",
	         "the reward baseline to separate successful and failed calls",
	         &heurdata->rewardbaseline, TRUE, DEFAULT_REWARDBASELINE, 0.0, 0.99, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/resetweights",
	         "should the bandit algorithms be reset when a new problem is read?",
	         &heurdata->resetweights, TRUE, DEFAULT_RESETWEIGHTS, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddStringParam(scip, "heuristics/" HEUR_NAME "/rewardfilename", "file name to store all rewards and the selection of the bandit",
	         &heurdata->rewardfilename, TRUE, DEFAULT_REWARDFILENAME, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/subsciprandseeds",
	         "should random seeds of sub-SCIPs be altered to increase diversification?",
	         &heurdata->subsciprandseeds, TRUE, DEFAULT_SUBSCIPRANDSEEDS, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/scalebyeffort",
	         "should the reward be scaled by the effort?",
	         &heurdata->scalebyeffort, TRUE, DEFAULT_SCALEBYEFFORT, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/copycuts",
	         "should cutting planes be copied to the sub-SCIP?",
	         &heurdata->copycuts, TRUE, DEFAULT_COPYCUTS, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/fixtol",
	         "tolerance by which the fixing rate may be missed without generic fixing",
	         &heurdata->fixtol, TRUE, DEFAULT_FIXTOL, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddRealParam(scip, "heuristics/" HEUR_NAME "/unfixtol",
	         "tolerance by which the fixing rate may be exceeded without generic unfixing",
	         &heurdata->unfixtol, TRUE, DEFAULT_UNFIXTOL, 0.0, 1.0, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/uselocalredcost",
	         "should local reduced costs be used for generic (un)fixing?",
	         &heurdata->uselocalredcost, TRUE, DEFAULT_USELOCALREDCOST, NULL, NULL) );
	
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/usepscost",
	         "should pseudo cost scores be used for variable priorization?",
	         &heurdata->usepscost, TRUE, DEFAULT_USEPSCOST, NULL, NULL) );
	   SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/" HEUR_NAME "/initduringroot",
	         "should the heuristic be executed multiple times during the root node?",
	         &heurdata->initduringroot, TRUE, DEFAULT_INITDURINGROOT, NULL, NULL) );
	
	   assert(SCIPfindTable(scip, TABLE_NAME_NEIGHBORHOOD) == NULL);
	   SCIP_CALL( SCIPincludeTable(scip, TABLE_NAME_NEIGHBORHOOD, TABLE_DESC_NEIGHBORHOOD, TRUE,
	         NULL, NULL, NULL, NULL, NULL, NULL, tableOutputNeighborhood,
	         NULL, TABLE_POSITION_NEIGHBORHOOD, TABLE_EARLIEST_STAGE_NEIGHBORHOOD) );
	
	   return SCIP_OKAY;
	}