1    	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2    	/*                                                                           */
3    	/*                  This file is part of the program and library             */
4    	/*         SCIP --- Solving Constraint Integer Programs                      */
5    	/*                                                                           */
6    	/*  Copyright (c) 2002-2023 Zuse Institute Berlin (ZIB)                      */
7    	/*                                                                           */
8    	/*  Licensed under the Apache License, Version 2.0 (the "License");          */
9    	/*  you may not use this file except in compliance with the License.         */
10   	/*  You may obtain a copy of the License at                                  */
11   	/*                                                                           */
12   	/*      http://www.apache.org/licenses/LICENSE-2.0                           */
13   	/*                                                                           */
14   	/*  Unless required by applicable law or agreed to in writing, software      */
15   	/*  distributed under the License is distributed on an "AS IS" BASIS,        */
16   	/*  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17   	/*  See the License for the specific language governing permissions and      */
18   	/*  limitations under the License.                                           */
19   	/*                                                                           */
20   	/*  You should have received a copy of the Apache-2.0 license                */
21   	/*  along with SCIP; see the file LICENSE. If not visit scipopt.org.         */
22   	/*                                                                           */
23   	/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24   	
25   	/**@file   sepa_lagromory.c
26   	 * @ingroup DEFPLUGINS_SEPA
27   	 * @brief  Lagromory separator
28   	 * @author Suresh Bolusani
29   	 * @author Mark Ruben Turner
30   	 * @author Mathieu Besançon
31   	 *
32   	 * This separator is based on the following article that discusses Lagromory separation using the relax-and-cut
33   	 * framework. Multiple enhancements have been implemented on top of the basic algorithm described in the article.
34   	 *
35   	 * Fischetti M. and Salvagnin D. (2011).@n
36   	 * A relax-and-cut framework for Gomory mixed-integer cuts.@n
37   	 * Mathematical Programming Computation, 3, 79-102.
38   	 *
39   	 * Consider the following linear relaxation at a node:
40   	 *
41   	 * \f[
42   	 *    \begin{array}{rrl}
43   	 *       \min & c^T x &\\
44   	 *            & x & \in P,
45   	 *    \end{array}
46   	 * \f]
47   	 *
48   	 * where \f$P\f$ is the feasible region of the relaxation. Let the following be the cuts generated so far in the current
49   	 * separation round.
50   	 *
51   	 * \f[
52   	 *    {\alpha^i}^T x \leq \alpha^i_0, i = 1, 2, \hdots, M
53   	 * \f]
54   	 *
55   	 * Then, the following is the Lagrangian dual problem considered in the relax-and-cut framework used in the separator.
56   	 *
57   	 * \f[
58   	 *    z_D := \max\limits_{u \geq 0} \left\{L(u) := \min \left\{c^T x + \sum\limits_{i = 1}^{M} \left(u_i
59   	 *    \left({\alpha^i}^T x - \alpha^i_0\right) \right) \mid x \in P\right\} \right\},
60   	 * \f]
61   	 *
62   	 * where \f$u\f$ are the Lagrangian multipliers (referred to as \a dualvector in this separator) used for penalizing the
63   	 * violation of the generated cuts, and \f$z_D\f$ is the optimal objective value (which is approximated via \a ubparam in this separator).
64   	 * Then, the following are the steps of the relax-and-cut algorithm implemented in this separator.
65   	 *
66   	 * \begin{itemize}
67   	 *    \item Generate an initial pool of cuts to build the initial Lagrangian dual problem.
68   	 *    \item Select initial values for Lagrangian multipliers \f$u^0\f$ (e.g., all zeroes vector).
69   	 *    \item In the outer main loop \f$i\f$ of the algorithm:
70   	 *       \begin{enumerate}
71   	 *          \item Solve the Lagrangian dual problem until certain termination criterion is met. This results in an inner
72   	 *          subgradient loop, whose iteration \f$j\f$ is described below.
73   	 *             \begin{enumerate}
74   	 *                \item Fix \f$u^j\f$, and solve the LP corresponding to the Lagrangian dual with fixed multipliers.
75   	 *                Gather its optimal simplex tableau and optimal objective value (i.e., the Lagrangian value)
76   	 *                \f$L(u^j)\f$.
77   	 *                \item Update \f$u^j\f$ to \f$u^{j+1}\f$ as follows.
78   	 *                   \f[
79   	 *                      u^{j+1} = \left(u^j + \lambda_j s^k\right)_+,
80   	 *                   \f]
81   	 *                   where \f$\lambda_j\f$ is the step length:
82   	 *                   \f[
83   	 *                      \lambda_j = \frac{\mu_j (UB - L(u^j))}{\|s^j\|^2_2},
84   	 *                   \f]
85   	 *                   where \f$mu_j\f$ is a factor (i.e., \a muparam) such that \f$0 < \mu_j \leq 2\f$, UB is \p ubparam,
86   	 *                   and \f$s^j\f$ is the subgradient vector defined as:
87   	 *                   \f[
88   	 *                      s^j_k = \left({\alpha^k}^T x - \alpha^k_0\right), k = 1, 2, \hdots, M.
89   	 *                   \f]
90   	 *                   The factor \f$mu_j\f$ is updated as below.
91   	 *                   \f[
92   	 *                      mu_j = \begin{cases}
93   	 *                               \p mubacktrackfactor * mu_j & \text{if } L(u^j) < bestLB - \delta\\
94   	 *                               \begin{cases}
95   	 *                                  \p muslab1factor * mu_j & \text{if } bestLB - avgLB < \p deltaslab1ub * delta\\
96   	 *                                  \p muslab2factor * mu_j & \text{if } \p deltaslab1ub * \delta \leq bestLB - avgLB < \p deltaslab2ub * delta\\
97   	 *                                  \p muslab3factor * mu_j & \text{otherwise}
98   	 *                               \end{cases} & \text{otherwise},
99   	 *                             \end{cases}
100  	 *                   \f]
101  	 *                   where \f$bestLB\f$ and \f$avgLB\f$ are best and average Lagrangian values found so far, and
102  	 *                   \f$\delta = UB - bestLB\f$.
103  	 *                \item Stabilize \f$u^{j+1}\f$ by projecting onto a norm ball followed by taking a convex combination
104  	 *                with a core vector of Lagrangian multipliers.
105  	 *                \item Generate GMI cuts based on the optimal simplex tableau.
106  	 *                \item Relax the newly generated cuts by penalizing and adding them to the objective function.
107  	 *                \item Go to the next iteration \f$j+1\f$.
108  	 *             \end{enumerate}
109  	 *          \item Gather all the generated cuts and build an LP by adding all these cuts to the node relaxation.
110  	 *          \item Solve this LP to obtain its optimal primal and dual solutions.
111  	 *          \item If this primal solution is MIP primal feasible, then add this solution to the solution pool, add all
112  	 *          the generated cuts to the cutpool or sepastore as needed, and exit the separator.
113  	 *          \item Otherwise, update the Lagrangian multipliers based on this optimal dual solution, and go to the next
114  	 *          iteration \f$i+1\f$.
115  	 *       \end{enumerate}
116  	 * \end{itemize}
117  	 *
118  	 * @todo store all LP sols in a data structure, and send them to fix-and-propagate at the end.
119  	 *
120  	 * @todo test heuristics such as feasibility pump with multiple input solutions.
121  	 *
122  	 * @todo find dual degenerate problems and test the separator on these problems.
123  	 *
124  	 * @todo identify instance classes where these cuts work better.
125  	 *
126  	 * @todo add termination criteria based on failed efforts.
127  	 *
128  	 * @todo for warm starting, if there are additional rows/cols, set their basis status to non-basic and then set WS info.
129  	 *
130  	 * @todo filter cuts using multiple explored LP solutions.
131  	 *
132  	 * @todo for bases on optimal face only, aggregate to get a new basis and separate it.
133  	 *
134  	 * @todo generate other separators in addition to GMI cuts (0-1/2)
135  	 *
136  	 * @todo: convert iters from int to SCIP_Longint
137  	 */
138  	
139  	/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
140  	
141  	#include "scip/heur_trysol.h"
142  	#include "scip/sepa_lagromory.h"
143  	
144  	#define SEPA_NAME                              "lagromory"
145  	#define SEPA_DESC                              "separator for Lagromory cuts for MIP relaxations"
146  	#define SEPA_PRIORITY                                -8000
147  	#define SEPA_FREQ                                       -1
148  	#define SEPA_MAXBOUNDDIST                              1.0
149  	#define SEPA_USESSUBSCIP                             FALSE /**< does the separator use a secondary SCIP instance? */
150  	#define SEPA_DELAY                                   FALSE /**< should separation method be delayed, if other separators found cuts? */
151  	
152  	/* generic parameters */
153  	#define DEFAULT_AWAY                                  0.01 /**< minimal integrality violation of a basis variable to try separation */
154  	#define DEFAULT_DELAYEDCUTS                          FALSE /**< should cuts be added to the delayed cut pool? */
155  	#define DEFAULT_SEPARATEROWS                          TRUE /**< separate rows with integral slack? */
156  	#define DEFAULT_SORTCUTOFFSOL                         TRUE /**< sort fractional integer columns based on fractionality? */
157  	#define DEFAULT_SIDETYPEBASIS                         TRUE /**< choose side types of row (lhs/rhs) based on basis information? */
158  	#define DEFAULT_DYNAMICCUTS                           TRUE /**< should generated cuts be removed from the LP if they are no longer tight? */
159  	#define DEFAULT_MAKEINTEGRAL                         FALSE /**< try to scale all cuts to integral coefficients? */
160  	#define DEFAULT_FORCECUTS                            FALSE /**< force cuts to be added to the LP? */
161  	#define DEFAULT_ALLOWLOCAL                           FALSE /**< should locally valid cuts be generated? */
162  	
163  	/* parameters related to the separator's termination check */
164  	#define DEFAULT_MAXROUNDSROOT                            1 /**< maximal number of separation rounds in the root node (-1: unlimited) */
165  	#define DEFAULT_MAXROUNDS                                1 /**< maximal number of separation rounds per node (-1: unlimited) */
166  	#define DEFAULT_DUALDEGENERACYRATETHRESHOLD            0.5 /**< minimum dual degeneracy rate for separator execution */
167  	#define DEFAULT_VARCONSRATIOTHRESHOLD                  1.0 /**< minimum variable-constraint ratio on optimal face for separator execution */
168  	#define DEFAULT_MINRESTART                               1 /**< minimum restart round for separator execution (0: from beginning of the
169  	                                                             instance solving, >= n with n >= 1: from restart round n) */
170  	#define DEFAULT_PERLPMAXCUTSROOT                        50 /**< maximal number of cuts separated per Lagromory LP in the root node */
171  	#define DEFAULT_PERLPMAXCUTS                            10 /**< maximal number of cuts separated per Lagromory LP in the non-root node */
172  	#define DEFAULT_PERROUNDLPITERLIMITFACTOR             -1.0 /**< factor w.r.t. root node LP iterations for maximal separating LP iterations
173  	                                                             per separation round (negative for no limit) */
174  	#define DEFAULT_ROOTLPITERLIMITFACTOR                 -1.0 /**< factor w.r.t. root node LP iterations for maximal separating LP iterations
175  	                                                             in the root node (negative for no limit) */
176  	#define DEFAULT_TOTALLPITERLIMITFACTOR                -1.0 /**< factor w.r.t. root node LP iterations for maximal separating LP iterations
177  	                                                             in the tree (negative for no limit) */
178  	#define DEFAULT_PERROUNDMAXLPITERS                   50000 /**< maximal number of separating LP iterations per separation round (-1: unlimited) */
179  	#define DEFAULT_PERROUNDCUTSFACTORROOT                 1.0 /**< factor w.r.t. number of integer columns for number of cuts separated per
180  	                                                             separation round in root node */
181  	#define DEFAULT_PERROUNDCUTSFACTOR                     0.5 /**< factor w.r.t. number of integer columns for number of cuts separated per
182  	                                                             separation round at a non-root node */
183  	#define DEFAULT_TOTALCUTSFACTOR                       50.0 /**< factor w.r.t. number of integer columns for total number of cuts separated */
184  	#define DEFAULT_MAXMAINITERS                             4 /**< maximal number of main loop iterations of the relax-and-cut algorithm */
185  	#define DEFAULT_MAXSUBGRADIENTITERS                      6 /**< maximal number of subgradient loop iterations of the relax-and-cut algorithm */
186  	
187  	/* parameters related to the relax-and-cut algorithm */
188  	#define DEFAULT_MUPARAMCONST                          TRUE /**< is the mu parameter (factor for step length) constant? */
189  	#define DEFAULT_MUPARAMINIT                           0.01 /**< initial value of the mu parameter (factor for step length) */
190  	#define DEFAULT_MUPARAMLB                              0.0 /**< lower bound for the mu parameter (factor for step length) */
191  	#define DEFAULT_MUPARAMUB                              2.0 /**< upper bound for the mu parameter (factor for step length) */
192  	#define DEFAULT_MUBACKTRACKFACTOR                      0.5 /**< factor of mu while backtracking the mu parameter (factor for step length) -
193  	                                                             see updateMuSteplengthParam() */
194  	#define DEFAULT_MUSLAB1FACTOR                         10.0 /**< factor of mu parameter (factor for step length) for larger increment - see
195  	                                                             updateMuSteplengthParam() */
196  	#define DEFAULT_MUSLAB2FACTOR                          2.0 /**< factor of mu parameter (factor for step length) for smaller increment - see
197  	                                                             updateMuSteplengthParam() */
198  	#define DEFAULT_MUSLAB3FACTOR                          0.5 /**< factor of mu parameter (factor for step length) for reduction - see
199  	                                                             updateMuSteplengthParam() */
200  	#define DEFAULT_DELTASLAB1UB                         0.001 /**< factor of delta deciding larger increment of mu parameter (factor for step
201  	                                                             length) - see updateMuSteplengthParam() */
202  	#define DEFAULT_DELTASLAB2UB                          0.01 /**< factor of delta deciding smaller increment of mu parameter (factor for step
203  	                                                             length) - see updateMuSteplengthParam() */
204  	#define DEFAULT_UBPARAMPOSFACTOR                       2.0 /**< factor for positive upper bound used as an estimate for the optimal
205  	                                                             Lagrangian dual value */
206  	#define DEFAULT_UBPARAMNEGFACTOR                       0.5 /**< factor for negative upper bound used as an estimate for the optimal
207  	                                                             Lagrangian dual value */
208  	#define DEFAULT_MAXLAGRANGIANVALSFORAVG                  2 /**< maximal number of iterations for rolling average of Lagrangian value */
209  	#define DEFAULT_MAXCONSECITERSFORMUUPDATE               10 /**< consecutive number of iterations used to determine if mu needs to be backtracked */
210  	#define DEFAULT_PERROOTLPITERFACTOR                    0.2 /**< factor w.r.t. root node LP iterations for iteration limit of each separating
211  	                                                             LP (negative for no limit) */
212  	#define DEFAULT_PERLPITERFACTOR                        0.1 /**< factor w.r.t. non-root node LP iterations for iteration limit of each
213  	                                                             separating LP (negative for no limit) */
214  	#define DEFAULT_CUTGENFREQ                               1 /**< frequency of subgradient iterations for generating cuts */
215  	#define DEFAULT_CUTADDFREQ                               1 /**< frequency of subgradient iterations for adding cuts to objective function */
216  	#define DEFAULT_CUTSFILTERFACTOR                       1.0 /**< fraction of generated cuts per explored basis to accept from separator */
217  	#define DEFAULT_OPTIMALFACEPRIORITY                      2 /**< priority of the optimal face for separator execution (0: low priority, 1:
218  	                                                             medium priority, 2: high priority) */
219  	#define DEFAULT_AGGREGATECUTS                         TRUE /**< aggregate all generated cuts using the Lagrangian multipliers? */
220  	/* parameters for stabilization of the Lagrangian multipliers */
221  	#define DEFAULT_PROJECTIONTYPE                           2 /**< the ball into which the Lagrangian multipliers are projected for
222  	                                                             stabilization (0: no projection, 1: L1-norm ball projection, 2: L2-norm ball
223  	                                                             projection, 3: L_inf-norm ball projection) */
224  	#define DEFAULT_STABILITYCENTERTYPE                      1 /**< type of stability center for taking weighted average of Lagrangian multipliers for
225  	                                                             stabilization (0: no weighted stabilization, 1: best Lagrangian multipliers) */
226  	#define DEFAULT_RADIUSINIT                             0.5 /**< initial radius of the ball used in stabilization of Lagrangian multipliers */
227  	#define DEFAULT_RADIUSMAX                             20.0 /**< maximum radius of the ball used in stabilization of Lagrangian multipliers */
228  	#define DEFAULT_RADIUSMIN                             1e-6 /**< minimum radius of the ball used in stabilization of Lagrangian multipliers */
229  	#define DEFAULT_CONST                                  2.0 /**< a constant for stablity center based stabilization of Lagrangian multipliers */
230  	#define DEFAULT_RADIUSUPDATEWEIGHT                    0.98 /**< multiplier to evaluate cut violation score used for updating ball radius */
231  	
232  	/* macros that are used directly */
233  	#define RANDSEED                                        42 /**< random seed */
234  	#define MAKECONTINTEGRAL                             FALSE /**< convert continuous variable to integral variables in SCIPmakeRowIntegral()? */
235  	#define POSTPROCESS                                   TRUE /**< apply postprocessing after MIR calculation? - see SCIPcalcMIR() */
236  	#define BOUNDSWITCH                                 0.9999 /**< threshold for bound switching - see SCIPcalcMIR() */
237  	#define USEVBDS                                       TRUE /**< use variable bounds? - see SCIPcalcMIR() */
238  	#define FIXINTEGRALRHS                               FALSE /**< try to generate an integral rhs? - see SCIPcalcMIR() */
239  	#define MAXAGGRLEN(ncols)               (0.1*(ncols)+1000) /**< maximal length of base inequality */
240  	
241  	/*
242  	 * Data structures
243  	 */
244  	
245  	/** separator data */
246  	struct SCIP_SepaData
247  	{
248  	   /* generic variables */
249  	   SCIP_Real             away;                        /**< minimal integrality violation of a basis variable to try separation */
250  	   SCIP_Bool             delayedcuts;                 /**< should cuts be added to the delayed cut pool? */
251  	   SCIP_Bool             separaterows;                /**< separate rows with integral slack? */
252  	   SCIP_Bool             sortcutoffsol;               /**< sort fractional integer columns based on fractionality? */
253  	   SCIP_Bool             sidetypebasis;               /**< choose side types of row (lhs/rhs) based on basis information? */
254  	   SCIP_Bool             dynamiccuts;                 /**< should generated cuts be removed from the LP if they are no longer tight? */
255  	   SCIP_Bool             makeintegral;                /**< try to scale all cuts to integral coefficients? */
256  	   SCIP_Bool             forcecuts;                   /**< force cuts to be added to the LP? */
257  	   SCIP_Bool             allowlocal;                  /**< should locally valid cuts be generated? */
258  	   SCIP_RANDNUMGEN*      randnumgen;                  /**< random number generator */
259  	   SCIP_HEUR*            heurtrysol;                  /**< a pointer to the trysol heuristic, if available */
260  	
261  	   /* used to define separating LPs */
262  	   SCIP_LPI*             lpiwithsoftcuts;             /**< pointer to the lpi interface of Lagrangian dual with fixed multipliers */
263  	   SCIP_LPI*             lpiwithhardcuts;             /**< pointer to the lpi interface of LP with all generated cuts */
264  	   int                   nrowsinhardcutslp;           /**< nrows of \a lpiwithhardcuts */
265  	   int                   nrunsinsoftcutslp;           /**< number of branch-and-bound runs on current instance */
266  	
267  	   /* used for termination checks */
268  	   SCIP_Longint          ncalls;                      /**< number of calls to the separator */
269  	   int                   maxroundsroot;               /**< maximal number of separation rounds in the root node (-1: unlimited) */
270  	   int                   maxrounds;                   /**< maximal number of separation rounds per node (-1: unlimited) */
271  	   SCIP_Real             dualdegeneracyratethreshold; /**< minimum dual degeneracy rate for separator execution */
272  	   SCIP_Real             varconsratiothreshold;       /**< minimum variable-constraint ratio on optimal face for separator execution */
273  	   int                   minrestart;                  /**< minimum restart round for separator execution (0: from beginning of the instance
274  	                                                        solving, >= n with n >= 1: from restart round n) */
275  	   int                   nmaxcutsperlproot;           /**< maximal number of cuts separated per Lagromory LP in the root node */
276  	   int                   nmaxcutsperlp;               /**< maximal number of cuts separated per Lagromory LP in the non-root node */
277  	   SCIP_Real             perroundlpiterlimitfactor;   /**< factor w.r.t. root node LP iterations for maximal separating LP iterations per
278  	                                                        separation round (negative for no limit) */
279  	   SCIP_Real             rootlpiterlimitfactor;       /**< factor w.r.t. root node LP iterations for maximal separating LP iterations in the
280  	                                                        root node (negative for no limit) */
281  	   SCIP_Real             totallpiterlimitfactor;      /**< factor w.r.t. root node LP iterations for maximal separating LP iterations in the
282  	                                                        tree (negative for no limit) */
283  	   int                   perroundnmaxlpiters;         /**< maximal number of separating LP iterations per separation round (-1: unlimited) */
284  	   SCIP_Real             perroundcutsfactorroot;      /**< factor w.r.t. number of integer columns for number of cuts separated per
285  	                                                        separation round in root node */
286  	   SCIP_Real             perroundcutsfactor;          /**< factor w.r.t. number of integer columns for number of cuts separated per
287  	                                                        separation round at a non-root node */
288  	   SCIP_Real             totalcutsfactor;             /**< factor w.r.t. number of integer columns for total number of cuts separated */
289  	   int                   nmaxmainiters;               /**< maximal number of main loop iterations of the relax-and-cut algorithm */
290  	   int                   nmaxsubgradientiters;        /**< maximal number of subgradient loop iterations of the relax-and-cut algorithm */
291  	   int                   nmaxperroundlpiters;         /**< maximal number of separating LP iterations per separation round */
292  	   int                   nmaxrootlpiters;             /**< maximal number of separating LP iterations in the root node */
293  	   int                   nrootlpiters;                /**< number of separating LP iterations in the root node */
294  	   int                   nmaxtotallpiters;            /**< maximal number of separating LP iterations in the tree */
295  	   int                   ntotallpiters;               /**< number of separating LP iterations in the tree */
296  	   int                   nmaxperroundcutsroot;        /**< maximal number of cuts separated per separation round in root node */
297  	   int                   nmaxperroundcuts;            /**< maximal number of cuts separated per separation round */
298  	   int                   nmaxtotalcuts;               /**< maximal number of cuts separated in the tree */
299  	   int                   ntotalcuts;                  /**< number of cuts separated in the tree */
300  	
301  	   /* used for the relax-and-cut algorithm */
302  	   SCIP_Bool             muparamconst;                /**< is the mu parameter (factor for step length) constant? */
303  	   SCIP_Real             muparaminit;                 /**< initial value of the mu parameter (factor for step length) */
304  	   SCIP_Real             muparamlb;                   /**< lower bound for the mu parameter (factor for step length) */
305  	   SCIP_Real             muparamub;                   /**< upper bound for the mu parameter (factor for step length) */
306  	   SCIP_Real             mubacktrackfactor;           /**< factor of mu while backtracking the mu parameter (factor for step length) - see
307  	                                                        updateMuSteplengthParam() */
308  	   SCIP_Real             muslab1factor;               /**< factor of mu parameter (factor for step length) for larger increment - see
309  	                                                        updateMuSteplengthParam() */
310  	   SCIP_Real             muslab2factor;               /**< factor of mu parameter (factor for step length) for smaller increment - see
311  	                                                        updateMuSteplengthParam() */
312  	   SCIP_Real             muslab3factor;               /**< factor of mu parameter (factor for step length) for reduction - see updateMuSteplengthParam() */
313  	   SCIP_Real             deltaslab1ub;                /**< factor of delta deciding larger increment of mu parameter (factor for step
314  	                                                        length) - see updateMuSteplengthParam() */
315  	   SCIP_Real             deltaslab2ub;                /**< factor of delta deciding smaller increment of mu parameter (factor for step
316  	                                                        length) - see updateMuSteplengthParam() */
317  	   SCIP_Real             ubparamposfactor;            /**< factor for positive upper bound used as an estimate for the optimal Lagrangian
318  	                                                        dual value */
319  	   SCIP_Real             ubparamnegfactor;            /**< factor for negative upper bound used as an estimate for the optimal Lagrangian
320  	                                                        dual value */
321  	   int                   nmaxlagrangianvalsforavg;    /**< maximal number of iterations for rolling average of Lagrangian value */
322  	   int                   nmaxconsecitersformuupdate;  /**< consecutive number of iterations used to determine if mu needs to be backtracked */
323  	   SCIP_Real             perrootlpiterfactor;         /**< factor w.r.t. root node LP iterations for iteration limit of each separating LP
324  	                                                        (negative for no limit) */
325  	   SCIP_Real             perlpiterfactor;             /**< factor w.r.t. non-root node LP iterations for iteration limit of each separating
326  	                                                        LP (negative for no limit) */
327  	   int                   cutgenfreq;                  /**< frequency of subgradient iterations for generating cuts */
328  	   int                   cutaddfreq;                  /**< frequency of subgradient iterations for adding cuts to objective function */
329  	   SCIP_Real             cutsfilterfactor;            /**< fraction of generated cuts per explored basis to accept from separator */
330  	   int                   optimalfacepriority;         /**< priority of the optimal face for separator execution (0: low priority, 1: medium
331  	                                                        priority, 2: high priority) */
332  	   SCIP_Bool             aggregatecuts;               /**< aggregate all generated cuts using the Lagrangian multipliers? */
333  	
334  	   /* for stabilization of Lagrangian multipliers */
335  	   int                   projectiontype;              /**< the ball into which the Lagrangian multipliers are projected for stabilization
336  	                                                        (0: no projection, 1: L1-norm ball projection, 2: L2-norm ball projection, 3:
337  	                                                        L_inf-norm ball projection) */
338  	   int                   stabilitycentertype;         /**< type of stability center for taking weighted average of Lagrangian multipliers for
339  	                                                        stabilization (0: no weighted stabilization, 1: best Lagrangian multipliers) */
340  	   SCIP_Real             radiusinit;                  /**< initial radius of the ball used in stabilization of Lagrangian multipliers */
341  	   SCIP_Real             radiusmax;                   /**< maximum radius of the ball used in stabilization of Lagrangian multipliers */
342  	   SCIP_Real             radiusmin;                   /**< minimum radius of the ball used in stabilization of Lagrangian multipliers */
343  	   SCIP_Real             constant;                    /**< a constant for stablity center based stabilization of Lagrangian multipliers */
344  	   SCIP_Real             radiusupdateweight;          /**< multiplier to evaluate cut violation score used for updating ball radius */
345  	};
346  	
347  	
348  	/*
349  	 * Local methods
350  	 */
351  	
352  	/** start the diving mode for solving LPs corresponding to the Lagrangian dual with fixed multipliers in the subgradient
353  	 * loop of the separator, and update some sepadata values */
354  	static
355  	SCIP_RETCODE createLPWithSoftCuts(
356  	   SCIP*                 scip,               /**< SCIP data structure */
357  	   SCIP_SEPADATA*        sepadata            /**< separator data structure */
358  	   )
359  	{
360  	   SCIP_ROW** rows;
361  	   SCIP_COL** cols;
362  	   int nruns;
363  	   int runnum;
364  	   int nrows;
365  	   int ncols;
366  	   unsigned int nintcols;
367  	   SCIP_Longint nrootlpiters;
368  	
369  	   assert(scip != NULL);
370  	   assert(sepadata != NULL);
371  	
372  	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
373  	   SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
374  	   runnum = SCIPgetNRuns(scip);  /* current run number of SCIP (starts from 1) indicating restarts */
375  	   nruns = sepadata->nrunsinsoftcutslp; /* previous run number of SCIP in which the diving LP was created */
376  	   nintcols = 0;
377  	
378  	   /* start diving mode so that the diving LP can be used for all subgradient iterations */
379  	   SCIP_CALL( SCIPstartDive(scip) );
380  	
381  	   /* store LPI pointer to be able to use LPI functions directly (e.g., setting time limit) */
382  	   SCIP_CALL( SCIPgetLPI(scip, &(sepadata->lpiwithsoftcuts)) );
383  	
384  	   /* if called for the first time in a restart (including the initial run), set certain sepadata values */
385  	   if( nruns != runnum )
386  	   {
387  	      /* get number of LP iterations of root node's first LP solving */
388  	      nrootlpiters = SCIPgetNRootFirstLPIterations(scip);
389  	
390  	      /* calculate maximum number of LP iterations allowed for all separation calls in the root node */
391  	      if( (sepadata->rootlpiterlimitfactor >= 0.0) && !SCIPisInfinity(scip, sepadata->rootlpiterlimitfactor) )
392  	      {
393  	         sepadata->nmaxrootlpiters = (int)(sepadata->rootlpiterlimitfactor * nrootlpiters);
394  	      }
395  	      else
396  	      {
397  	         sepadata->nmaxrootlpiters = -1; /* no finite limit */
398  	      }
399  	
400  	      /* calculate maximum number of LP iterations allowed for all separation calls in the entire tree */
401  	      if( (sepadata->totallpiterlimitfactor >= 0.0) && !SCIPisInfinity(scip, sepadata->totallpiterlimitfactor) )
402  	      {
403  	         sepadata->nmaxtotallpiters = (int)(sepadata->totallpiterlimitfactor * nrootlpiters);
404  	      }
405  	      else
406  	      {
407  	         sepadata->nmaxtotallpiters = -1; /* no finite limit */
408  	      }
409  	
410  	      /* calculate maximum number of LP iterations allowed per separation call */
411  	      if( (sepadata->perroundlpiterlimitfactor >= 0.0) && !SCIPisInfinity(scip, sepadata->perroundlpiterlimitfactor) )
412  	      {
413  	         sepadata->nmaxperroundlpiters = (int)(sepadata->perroundlpiterlimitfactor * nrootlpiters);
414  	      }
415  	      else
416  	      {
417  	         sepadata->nmaxperroundlpiters = -1; /* no finite limit */
418  	      }
419  	
420  	      /* update maximum number of LP iterations allowed per separation call using absolute limits */
421  	      if( sepadata->perroundnmaxlpiters > 0 )
422  	      {
423  	         sepadata->nmaxperroundlpiters = ((sepadata->nmaxperroundlpiters >= 0) ? MIN(sepadata->nmaxperroundlpiters,
424  	                  sepadata->perroundnmaxlpiters) : sepadata->perroundnmaxlpiters);
425  	      }
426  	
427  	      /* set maximum number of cuts allowed to generate per round in root and non-root nodes as well as the total tree */
428  	      for( int i = 0; i < ncols; ++i )
429  	      {
430  	         nintcols += SCIPcolIsIntegral(cols[i]);
431  	      }
432  	      sepadata->nmaxperroundcutsroot = (int)(sepadata->perroundcutsfactorroot * nintcols);
433  	      sepadata->nmaxperroundcuts = (int)(sepadata->perroundcutsfactor * nintcols);
434  	
435  	      if( sepadata->ncalls == 0 )
436  	      {
437  	         sepadata->nmaxtotalcuts = (int)(sepadata->totalcutsfactor * nintcols);
438  	         sepadata->ntotalcuts = 0;
439  	      }
440  	
441  	      /* update the run number of solving to represent the restart number in which the above limits were set */
442  	      sepadata->nrunsinsoftcutslp = runnum;
443  	   }
444  	
445  	   return SCIP_OKAY;
446  	}
447  	
448  	/** end the diving mode that was used for solving LPs corresponding to the Lagrangian dual with fixed multipliers */
449  	static
450  	SCIP_RETCODE deleteLPWithSoftCuts(
451  	   SCIP*                 scip,               /**< SCIP data structure */
452  	   SCIP_SEPADATA*        sepadata            /**< separator data structure */
453  	   )
454  	{
455  	   assert(scip != NULL);
456  	   assert(sepadata != NULL);
457  	
458  	   SCIP_CALL( SCIPendDive(scip) );
459  	
460  	   return SCIP_OKAY;
461  	}
462  	
463  	/** set up LP interface to solve LPs in the (outer) main loop of the relax-and-cut algorithm; these LPs are built by
464  	 * adding all the generated cuts to the node relaxation */
465  	/* @todo add lpi iters to global statistics */
466  	static
467  	SCIP_RETCODE createLPWithHardCuts(
468  	   SCIP*                 scip,               /**< SCIP data structure */
469  	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
470  	   SCIP_ROW**            cuts,               /**< generated cuts to be added to the LP */
471  	   int                   ncuts               /**< number of generated cuts to be added to the LP */
472  	   )
473  	{
474  	   SCIP_ROW** rows;
475  	   SCIP_COL** cols;
476  	   SCIP_Real* collb;
477  	   SCIP_Real* colub;
478  	   SCIP_Real* colobj;
479  	   SCIP_Real* rowlhs;
480  	   SCIP_Real* rowrhs;
481  	   SCIP_Real rowconst;
482  	   SCIP_Real* rowvals;
483  	   SCIP_Real* rowval;
484  	   SCIP_COL** rowcols;
485  	   SCIP_LPI* wslpi;
486  	   SCIP_LPISTATE* lpistate;
487  	   BMS_BLKMEM* blkmem;
488  	   SCIP_Real pinf;
489  	   SCIP_Real ninf;
490  	   int nrows;
491  	   int ncols;
492  	   int nrownonz;
493  	   int collppos;
494  	   int* rowcolinds;
495  	   int* rowbegs;
496  	
497  	   assert(scip != NULL);
498  	   assert(sepadata != NULL);
499  	
500  	   SCIP_LPI** lpi = &(sepadata->lpiwithhardcuts);
501  	   blkmem = SCIPblkmem(scip);
502  	
503  	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
504  	   SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
505  	
506  	   /* if this function is called for the first time in this separation round, then create an LPI and add cols & rows */
507  	   if( ncuts == 0 )
508  	   {
509  	      if( *lpi != NULL )
510  	      {
511  	         SCIP_CALL( SCIPlpiFree(lpi) );
512  	         *lpi = NULL;
513  	      }
514  	      assert(*lpi == NULL);
515  	
516  	      /* create an LPI with appropriate objective sense */
517  	      if( SCIPgetObjsense(scip) == SCIP_OBJSENSE_MAXIMIZE )
518  	      {
519  	         SCIP_CALL( SCIPlpiCreate(lpi, SCIPgetMessagehdlr(scip), "node LP with generated cuts", SCIP_OBJSEN_MAXIMIZE) );
520  	      }
521  	      else
522  	      {
523  	         SCIP_CALL( SCIPlpiCreate(lpi, SCIPgetMessagehdlr(scip), "node LP with generated cuts", SCIP_OBJSEN_MINIMIZE) );
524  	      }
525  	
526  	      /* add cols to the LP interface */
527  	      SCIP_CALL( SCIPallocBufferArray(scip, &colobj, ncols) );
528  	      SCIP_CALL( SCIPallocBufferArray(scip, &collb, ncols) );
529  	      SCIP_CALL( SCIPallocBufferArray(scip, &colub, ncols) );
530  	      /* gather required column information */
531  	      for( int i = 0; i < ncols; ++i )
532  	      {
533  	         colobj[i] = SCIPcolGetObj(cols[i]);
534  	         collb[i] = SCIPcolGetLb(cols[i]);
535  	         colub[i] = SCIPcolGetUb(cols[i]);
536  	      }
537  	      /* add cols */
538  	      SCIP_CALL( SCIPlpiAddCols(*lpi, ncols, colobj, collb, colub, NULL, 0, NULL, NULL, NULL) );
539  	      SCIPfreeBufferArray(scip, &colub);
540  	      SCIPfreeBufferArray(scip, &collb);
541  	      SCIPfreeBufferArray(scip, &colobj);
542  	
543  	      /* add rows to the LP interface */
544  	      /* find number of nonzeroes in rows */
545  	      nrownonz = 0;
546  	      for( int i = 0; i < nrows; ++i )
547  	      {
548  	         assert(!(SCIPisInfinity(scip, -SCIProwGetLhs(rows[i])) && SCIPisInfinity(scip, SCIProwGetRhs(rows[i]))));
549  	         nrownonz += SCIProwGetNLPNonz(rows[i]);
550  	      }
551  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowcolinds, nrownonz) );
552  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowvals, nrownonz) );
553  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowbegs, nrows + 1) );
554  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowlhs, nrows) );
555  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowrhs, nrows) );
556  	      /* gather required row information */
557  	      rowbegs[0] = 0;
558  	      pinf = SCIPlpiInfinity(*lpi);
559  	      ninf = -SCIPlpiInfinity(*lpi);
560  	      for( int i = 0; i < nrows; ++i )
561  	      {
562  	         nrownonz = SCIProwGetNLPNonz(rows[i]);
563  	         assert(nrownonz <= ncols);
564  	         rowval = SCIProwGetVals(rows[i]);
565  	         rowcols = SCIProwGetCols(rows[i]);
566  	
567  	         rowbegs[i + 1] = rowbegs[i] + nrownonz;
568  	         rowconst = SCIProwGetConstant(rows[i]);
569  	         rowlhs[i] = SCIPisInfinity(scip, -SCIProwGetLhs(rows[i])) ? ninf : SCIProwGetLhs(rows[i]) - rowconst;
570  	         rowrhs[i] = SCIPisInfinity(scip, SCIProwGetRhs(rows[i])) ? pinf : SCIProwGetRhs(rows[i]) - rowconst;
571  	
572  	         for( int j = 0; j < nrownonz; ++j )
573  	         {
574  	            collppos = SCIPcolGetLPPos(rowcols[j]);
575  	            assert(collppos >= 0);
576  	            assert(collppos <= ncols);
577  	
578  	            rowcolinds[rowbegs[i] + j] = collppos;
579  	            rowvals[rowbegs[i] + j] = rowval[j];
580  	         }
581  	      }
582  	      /* add rows */
583  	      SCIP_CALL( SCIPlpiAddRows(*lpi, nrows, rowlhs, rowrhs, NULL, rowbegs[nrows], rowbegs, rowcolinds, rowvals) );
584  	
585  	      /* get warm starting info */
586  	      SCIP_CALL( SCIPgetLPI(scip, &wslpi) );
587  	      SCIP_CALL( SCIPlpiGetState(wslpi, blkmem, &lpistate) );
588  	   }
589  	   /* if there are any cuts, then add the cuts that were not added earlier to the LPI */
590  	   else
591  	   {
592  	      assert(nrows + ncuts >= sepadata->nrowsinhardcutslp);
593  	
594  	      /* get warm starting info */
595  	      wslpi = *lpi;
596  	      SCIP_CALL( SCIPlpiGetState(wslpi, blkmem, &lpistate) );
597  	
598  	      /* find number of nonzeros in cuts and allocate memory */
599  	      nrownonz = 0;
600  	      pinf = SCIPlpiInfinity(*lpi);
601  	      ninf = -SCIPlpiInfinity(*lpi);
602  	      for( int i = sepadata->nrowsinhardcutslp - nrows; i < ncuts; ++i )
603  	      {
604  	         assert(!(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])) && SCIPisInfinity(scip, SCIProwGetRhs(cuts[i]))));
605  	         assert(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])));
606  	         assert(!SCIPisInfinity(scip, SCIProwGetRhs(cuts[i])));
607  	         nrownonz += SCIProwGetNNonz(cuts[i]);
608  	      }
609  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowcolinds, nrownonz) );
610  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowvals, nrownonz) );
611  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowbegs, (ncuts - sepadata->nrowsinhardcutslp + nrows) + 1) );
612  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowlhs, (ncuts - sepadata->nrowsinhardcutslp + nrows)) );
613  	      SCIP_CALL( SCIPallocBufferArray(scip, &rowrhs, (ncuts - sepadata->nrowsinhardcutslp + nrows)) );
614  	
615  	      /* gather required cut information */
616  	      rowbegs[0] = 0;
617  	      for( int i = sepadata->nrowsinhardcutslp - nrows; i < ncuts; ++i )
618  	      {
619  	         nrownonz = SCIProwGetNNonz(cuts[i]);
620  	         assert(nrownonz <= ncols);
621  	         rowval = SCIProwGetVals(cuts[i]);
622  	         rowcols = SCIProwGetCols(cuts[i]);
623  	
624  	         rowbegs[i - sepadata->nrowsinhardcutslp + nrows + 1] = rowbegs[i - sepadata->nrowsinhardcutslp + nrows] +
625  	            nrownonz;
626  	         rowconst = SCIProwGetConstant(cuts[i]);
627  	         rowlhs[i - sepadata->nrowsinhardcutslp + nrows] = SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])) ? ninf :
628  	            SCIProwGetLhs(cuts[i]) - rowconst;
629  	         rowrhs[i - sepadata->nrowsinhardcutslp + nrows] = SCIPisInfinity(scip, SCIProwGetRhs(cuts[i])) ? pinf :
630  	            SCIProwGetRhs(cuts[i]) - rowconst;
631  	
632  	         for( int j = 0; j < nrownonz; ++j )
633  	         {
634  	            collppos = SCIPcolGetLPPos(rowcols[j]);
635  	            assert(collppos >= 0);
636  	            assert(collppos <= ncols);
637  	
638  	            rowcolinds[rowbegs[i - sepadata->nrowsinhardcutslp + nrows] + j] = collppos;
639  	            rowvals[rowbegs[i - sepadata->nrowsinhardcutslp + nrows] + j] = rowval[j];
640  	         }
641  	      }
642  	
643  	      /* add cuts */
644  	      SCIP_CALL( SCIPlpiAddRows(*lpi, (ncuts - sepadata->nrowsinhardcutslp + nrows), rowlhs, rowrhs, NULL,
645  	               rowbegs[(ncuts - sepadata->nrowsinhardcutslp + nrows)], rowbegs, rowcolinds, rowvals) );
646  	   }
647  	
648  	   /* set warm starting basis */
649  	   SCIP_CALL( SCIPlpiSetState(*lpi, blkmem, lpistate) );
650  	
651  	   /* reset remaining sepadata values */
652  	   sepadata->nrowsinhardcutslp = nrows + ncuts;
653  	
654  	   /* free memory */
655  	   SCIP_CALL( SCIPlpiFreeState(*lpi, blkmem, &lpistate) );
656  	   SCIPfreeBufferArray(scip, &rowrhs);
657  	   SCIPfreeBufferArray(scip, &rowlhs);
658  	   SCIPfreeBufferArray(scip, &rowbegs);
659  	   SCIPfreeBufferArray(scip, &rowvals);
660  	   SCIPfreeBufferArray(scip, &rowcolinds);
661  	
662  	   return SCIP_OKAY;
663  	}
664  	
665  	/** free separator data */
666  	static
667  	SCIP_RETCODE sepadataFree(
668  	   SCIP*                 scip,               /**< SCIP data structure */
669  	   SCIP_SEPADATA**       sepadata            /**< separator data structure */
670  	   )
671  	{
672  	   assert(scip != NULL);
673  	   assert(sepadata != NULL);
674  	   assert(*sepadata != NULL);
675  	
676  	   if( (*sepadata)->lpiwithhardcuts != NULL )
677  	   {
678  	      SCIP_CALL( SCIPlpiFree(&((*sepadata)->lpiwithhardcuts)) );
679  	   }
680  	
681  	   (*sepadata)->nrowsinhardcutslp = 0;
682  	   (*sepadata)->nrunsinsoftcutslp = 0;
683  	   (*sepadata)->ncalls = 0;
684  	   (*sepadata)->nmaxperroundlpiters = 0;
685  	   (*sepadata)->nmaxrootlpiters = 0;
686  	   (*sepadata)->nrootlpiters = 0;
687  	   (*sepadata)->nmaxtotallpiters = 0;
688  	   (*sepadata)->ntotallpiters = 0;
689  	   (*sepadata)->nmaxperroundcutsroot = 0;
690  	   (*sepadata)->nmaxperroundcuts = 0;
691  	   (*sepadata)->nmaxtotalcuts = 0;
692  	   (*sepadata)->ntotalcuts = 0;
693  	
694  	   SCIPfreeBlockMemory(scip, sepadata);
695  	
696  	   return SCIP_OKAY;
697  	}
698  	
699  	/** update mu parameter which is used as a factor in the step length calculation; refer to the top of the file for a
700  	 * description of the formula.
701  	 */
702  	/* @todo some adaptive strategy like constant after certain changes? */
703  	static
704  	SCIP_RETCODE updateMuSteplengthParam(
705  	   SCIP*                 scip,               /**< SCIP data structure */
706  	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
707  	   int                   subgradientiternum, /**< subgradient iteration number */
708  	   SCIP_Real             ubparam,            /**< estimate of the optimal Lagrangian dual value */
709  	   SCIP_Real*            lagrangianvals,     /**< vector of Lagrangian values found so far */
710  	   SCIP_Real             bestlagrangianval,  /**< best Lagrangian value found so far */
711  	   SCIP_Real             avglagrangianval,   /**< rolling average of the Lagrangian values found so far */
712  	   SCIP_Real*            muparam,            /**< mu parameter to be updated */
713  	   SCIP_Bool*            backtrack           /**< whether mu parameter has been backtracked */
714  	   )
715  	{
716  	   SCIP_Real delta;
717  	   SCIP_Real deltaslab1ub;
718  	   SCIP_Real deltaslab2ub;
719  	   SCIP_Real muslab1factor;
720  	   SCIP_Real muslab2factor;
721  	   SCIP_Real muslab3factor;
722  	   int maxiters;
723  	   int i;
724  	
725  	   *backtrack = FALSE;
726  	
727  	   /* update the mu parameter only if it is not set to be a constant value */
728  	   if( !sepadata->muparamconst )
729  	   {
730  	      delta = ubparam - bestlagrangianval;
731  	      deltaslab1ub = MIN(sepadata->deltaslab1ub, sepadata->deltaslab2ub);
732  	      deltaslab2ub = MAX(sepadata->deltaslab1ub, sepadata->deltaslab2ub);
733  	      /* ensure that the ordering of different user input parameters is as expected */
734  	      if( SCIPisPositive(scip, sepadata->muslab1factor - sepadata->muslab2factor) )
735  	      {
736  	         if( SCIPisPositive(scip, sepadata->muslab2factor - sepadata->muslab3factor) )
737  	         {
738  	            muslab1factor = sepadata->muslab1factor;
739  	            muslab2factor = sepadata->muslab2factor;
740  	            muslab3factor = sepadata->muslab3factor;
741  	         }
742  	         else
743  	         {
744  	            if( SCIPisPositive(scip, sepadata->muslab1factor - sepadata->muslab3factor) )
745  	            {
746  	               muslab1factor = sepadata->muslab1factor;
747  	               muslab2factor = sepadata->muslab3factor;
748  	               muslab3factor = sepadata->muslab2factor;
749  	            }
750  	            else
751  	            {
752  	               muslab1factor = sepadata->muslab3factor;
753  	               muslab2factor = sepadata->muslab1factor;
754  	               muslab3factor = sepadata->muslab2factor;
755  	            }
756  	         }
757  	      }
758  	      else
759  	      {
760  	         if( SCIPisPositive(scip, sepadata->muslab1factor - sepadata->muslab3factor) )
761  	         {
762  	            muslab1factor = sepadata->muslab2factor;
763  	            muslab2factor = sepadata->muslab1factor;
764  	            muslab3factor = sepadata->muslab3factor;
765  	         }
766  	         else
767  	         {
768  	            if( SCIPisPositive(scip, sepadata->muslab2factor - sepadata->muslab3factor) )
769  	            {
770  	               muslab1factor = sepadata->muslab2factor;
771  	               muslab2factor = sepadata->muslab3factor;
772  	               muslab3factor = sepadata->muslab1factor;
773  	            }
774  	            else
775  	            {
776  	               muslab1factor = sepadata->muslab3factor;
777  	               muslab2factor = sepadata->muslab2factor;
778  	               muslab3factor = sepadata->muslab1factor;
779  	            }
780  	         }
781  	      }
782  	
783  	      maxiters = MIN(sepadata->nmaxconsecitersformuupdate, sepadata->nmaxlagrangianvalsforavg);
784  	      i = -1;
785  	
786  	      /* if certain number of iterations are done, then check for a possibility of backtracking and apply accordingly */
787  	      if( subgradientiternum >= maxiters )
788  	      {
789  	         for( i = subgradientiternum - maxiters; i < subgradientiternum; i++ )
790  	         {
791  	            if( SCIPisGE(scip, lagrangianvals[i], bestlagrangianval - delta) )
792  	               break;
793  	         }
794  	
795  	         if( i == subgradientiternum )
796  	         {
797  	            *muparam *= sepadata->mubacktrackfactor;
798  	            *backtrack = TRUE;
799  	         }
800  	      }
801  	
802  	      /* update mu parameter based on the different between best and average Lagrangian values */
803  	      if( (subgradientiternum < maxiters) || (i >= 0 && i < subgradientiternum) )
804  	      {
805  	         if( bestlagrangianval - avglagrangianval < deltaslab1ub * delta )
806  	            *muparam *= muslab1factor;
807  	         else if( bestlagrangianval - avglagrangianval < deltaslab2ub * delta )
808  	            *muparam *= muslab2factor;
809  	         else
810  	            *muparam *= muslab3factor;
811  	      }
812  	
813  	      /* reset the mu parameter to within its bounds */
814  	      *muparam = MAX(*muparam, sepadata->muparamlb);
815  	      *muparam = MIN(*muparam, sepadata->muparamub);
816  	   }
817  	
818  	   return SCIP_OKAY;
819  	}
820  	
821  	/** update subgradient, i.e., residuals of generated cuts */
822  	/* @note: assumed that \f$i^{th}\f$ cut is of the form \f${\alpha^i}^T x \leq {\alpha^i_0}\f$ */
823  	static
824  	void updateSubgradient(
825  	   SCIP*                 scip,                     /**< SCIP data structure */
826  	   SCIP_SOL*             sol,                      /**< LP solution used in updating subgradient vector */
827  	   SCIP_ROW**            cuts,                     /**< cuts generated so far */
828  	   int                   ncuts,                    /**< number of cuts generated so far */
829  	   SCIP_Real*            subgradient,              /**< vector of subgradients to be updated */
830  	   SCIP_Real*            dualvector,               /**< Lagrangian multipliers */
831  	   SCIP_Bool*            subgradientzero,          /**< whether the subgradient vector is all zero */
832  	   int*                  ncutviols,                /**< number of violations of generated cuts */
833  	   SCIP_Real*            maxcutviol,               /**< maximum violation of generated cuts */
834  	   int*                  nnzsubgradientdualprod,   /**< number of nonzero products of subgradient vector and Lagrangian multipliers (i.e.,
835  	                                                     number of complementarity slackness violations) */
836  	   SCIP_Real*            maxnzsubgradientdualprod  /**< maximum value of nonzero products of subgradient vector and Lagrangian multipliers
837  	                                                     (i.e., maximum value of complementarity slackness violations) */
838  	   )
839  	{
840  	   int nzerosubgradient;
841  	   SCIP_Real term;
842  	
843  	   assert(subgradientzero != NULL);
844  	   assert(ncutviols != NULL);
845  	   assert(maxcutviol != NULL);
846  	   assert(nnzsubgradientdualprod != NULL);
847  	   assert(maxnzsubgradientdualprod != NULL);
848  	
849  	   *ncutviols = 0;
850  	   *maxcutviol = 0.0;
851  	   *nnzsubgradientdualprod = 0;
852  	   *maxnzsubgradientdualprod = 0.0;
853  	   nzerosubgradient = 0;
854  	   *subgradientzero = FALSE;
855  	
856  	   /* for each cut, calculate the residual along with various violation metrics */
857  	   for( int i = 0; i < ncuts; i++ )
858  	   {
859  	      assert(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])));
860  	      assert(!SCIPisInfinity(scip, SCIProwGetRhs(cuts[i])));
861  	      subgradient[i] = SCIPgetRowSolActivity(scip, cuts[i], sol) + SCIProwGetConstant(cuts[i]) - SCIProwGetRhs(cuts[i]);
862  	      if( SCIPisFeasZero(scip, subgradient[i]) )
863  	      {
864  	         subgradient[i] = 0.0;
865  	         nzerosubgradient++;
866  	      }
867  	      else
868  	      {
869  	         /* check for cut violation */
870  	         if( SCIPisFeasPositive(scip, subgradient[i]) )
871  	         {
872  	            (*ncutviols)++;
873  	            *maxcutviol = MAX(*maxcutviol, subgradient[i]);
874  	         }
875  	
876  	         /* check for violation of complementarity slackness associated with the cut */
877  	         if( !SCIPisZero(scip, subgradient[i] * dualvector[i]) )
878  	         {
879  	            (*nnzsubgradientdualprod)++;
880  	            term = REALABS(subgradient[i] * dualvector[i]);
881  	            *maxnzsubgradientdualprod = MAX(*maxnzsubgradientdualprod, term);
882  	         }
883  	      }
884  	   }
885  	
886  	   /* indicator for all zero subgradient vector */
887  	   if( nzerosubgradient == ncuts )
888  	   {
889  	      *subgradientzero = TRUE;
890  	   }
891  	}
892  	
893  	/** update Lagrangian value, i.e., optimal value of the Lagrangian dual with fixed multipliers */
894  	static
895  	void updateLagrangianValue(
896  	   SCIP*                 scip,               /**< SCIP data structure */
897  	   SCIP_Real             objval,             /**< objective value of the Lagrangian dual with fixed multipliers */
898  	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers */
899  	   SCIP_ROW**            cuts,               /**< cuts generated so far */
900  	   int                   ncuts,              /**< number of cuts generated so far */
901  	   SCIP_Real*            lagrangianval       /**< Lagrangian value to be updated */
902  	   )
903  	{
904  	   *lagrangianval = objval;
905  	
906  	   for( int i = 0; i < ncuts; i++ )
907  	   {
908  	      assert(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])));
909  	      assert(!SCIPisInfinity(scip, SCIProwGetRhs(cuts[i])));
910  	      *lagrangianval += dualvector[i] * (SCIProwGetConstant(cuts[i]) - SCIProwGetRhs(cuts[i]));
911  	   }
912  	}
913  	
914  	/** update step length based on various input arguments; refer to the top of the file for an expression */
915  	static
916  	void updateStepLength(
917  	   SCIP*                 scip,               /**< SCIP data structure */
918  	   SCIP_Real             muparam,            /**< mu parameter used as a factor for step length */
919  	   SCIP_Real             ubparam,            /**< estimate of the optimal Lagrangian dual value */
920  	   SCIP_Real             lagrangianval,      /**< Lagrangian value */
921  	   SCIP_Real*            subgradient,        /**< subgradient vector */
922  	   int                   ncuts,              /**< number of cuts generated so far */
923  	   SCIP_Real*            steplength          /**< step length to be updated */
924  	   )
925  	{
926  	   SCIP_Real normsquared = 0.0;
927  	
928  	   for( int i = 0; i < ncuts; i++ )
929  	   {
930  	      normsquared += SQR(subgradient[i]);
931  	   }
932  	
933  	   if( !SCIPisFeasZero(scip, normsquared) )
934  	   {
935  	      *steplength = (muparam * (ubparam - lagrangianval))/(normsquared); /*lint !e795*/
936  	   }
937  	}
938  	
939  	/** update the ball radius (based on various violation metrics) that is used for stabilization of Lagrangian multipliers */
940  	static
941  	void updateBallRadius(
942  	   SCIP*                 scip,               /**< SCIP data structure */
943  	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
944  	   SCIP_Real             maxviolscore,       /**< weighted average of maximum value of generated cut violations and maximum value of
945  	                                               complementarity slackness violations, in the current iteration */
946  	   SCIP_Real             maxviolscoreold,    /**< weighted average of maximum value of generated cut violations and maximum value of
947  	                                               complementarity slackness violations, in the previous iteration */
948  	   SCIP_Real             nviolscore,         /**< weighted average of number of generated cut violations and number of complementarity
949  	                                               slackness violations, in the current iteration */
950  	   SCIP_Real             nviolscoreold,      /**< weighted average of number of generated cut violations and number of complementarity
951  	                                               slackness violations, in the previous iteration */
952  	   int                   nlpiters,           /**< number of LP iterations taken for solving the Lagrangian dual with fixed multipliers in
953  	                                               current iteration */
954  	   SCIP_Real*            ballradius          /**< norm ball radius to be updated */
955  	   )
956  	{
957  	   SCIP_Bool maxviolscoreimproved;
958  	   SCIP_Bool nviolscoreimproved;
959  	
960  	   assert(ballradius != NULL);
961  	
962  	   maxviolscoreimproved = !SCIPisNegative(scip, maxviolscoreold - maxviolscore);
963  	   nviolscoreimproved = !SCIPisNegative(scip, nviolscoreold - nviolscore);
964  	
965  	   if( maxviolscoreimproved && nviolscoreimproved )
966  	   {
967  	      /* both the maximum violation and number of violations scores have become better, so, increase the radius */
968  	      if( sepadata->optimalfacepriority <= 1 )
969  	      {
970  	         *ballradius *= 2.0;
971  	         *ballradius = MIN(*ballradius, sepadata->radiusmax);
972  	      }
973  	      else
974  	      {
975  	         *ballradius *= 1.5;
976  	         *ballradius = MIN(*ballradius, sepadata->radiusmax/2.0);
977  	      }
978  	   }
979  	   else if( !maxviolscoreimproved && !nviolscoreimproved )
980  	   {
981  	      /* both the maximum violation and number of violations scores have become worse, so, decrease the radius */
982  	      *ballradius *= 0.5;
983  	      *ballradius = MAX(*ballradius, sepadata->radiusmin);
984  	   }
985  	   else if( nlpiters == 0 )
986  	   {
987  	      /* only one among the maximum violation and number of violations scores has become better, and the LP basis did
988  	       * not change (i.e., nlpters = 0), so, increase the radius slightly */
989  	      if( sepadata->optimalfacepriority <= 1 )
990  	      {
991  	         *ballradius *= 1.5;
992  	         *ballradius = MIN(*ballradius, sepadata->radiusmax);
993  	      }
994  	      else
995  	      {
996  	         *ballradius *= 1.2;
997  	         *ballradius = MIN(*ballradius, sepadata->radiusmax/2.0);
998  	      }
999  	   }
1000 	}
1001 	
1002 	/** projection of Lagrangian multipliers onto L1-norm ball. This algorithm is based on the following article.
1003 	 *
1004 	 * Condat L. (2016).@n
1005 	 * Fast projection onto the simplex and the \f$l_1\f$ ball.@n
1006 	 * Mathematical Programming, 158, 1-2, 575-585.
1007 	 *
1008 	 */
1009 	static
1010 	SCIP_RETCODE l1BallProjection(
1011 	   SCIP*                 scip,               /**< SCIP data structure */
1012 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers to be projected onto L1-norm vall */
1013 	   int                   dualvectorlen,      /**< length of the Lagrangian multipliers vector */
1014 	   SCIP_Real             radius              /**< radius of the L1-norm ball */
1015 	   )
1016 	{
1017 	   SCIP_Real* temp1vals;
1018 	   SCIP_Real* temp2vals;
1019 	   SCIP_Real pivotparam;
1020 	   SCIP_Real val;
1021 	   SCIP_Real term;
1022 	   SCIP_Bool temp1changed;
1023 	   int ntemp1removed;
1024 	   int* temp1inds;
1025 	   int* temp2inds;
1026 	   int temp1len;
1027 	   int temp2len;
1028 	
1029 	   assert(!SCIPisNegative(scip, radius));
1030 	   val = REALABS(dualvector[0]);
1031 	   /* calculate the L1-norm of the Lagrangian multipliers */
(1) Event cond_false: Condition "i < dualvectorlen", taking false branch.
1032 	   for( int i = 1; i < dualvectorlen; i++ )
1033 	   {
1034 	      val += REALABS(dualvector[i]);
(2) Event loop_end: Reached end of loop.
1035 	   }
1036 	
1037 	   /* if the vector of Lagrangian multipliers lies outside the L1-norm ball, then do the projection */
(3) Event cond_true: Condition "val - radius > scip->set->num_epsilon", taking true branch.
1038 	   if( SCIPisGT(scip, val, radius) )
1039 	   {
(4) Event cond_false: Condition "(temp1vals = BMSallocBufferMemoryArray_call(SCIPcleanbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 8UL /* sizeof (*temp1vals) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1040)) == NULL", taking false branch.
(5) Event cond_false: Condition "(_restat_ = (((temp1vals = BMSallocBufferMemoryArray_call(SCIPcleanbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 8UL /* sizeof (*temp1vals) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1040)) == NULL) ? SCIP_NOMEMORY : SCIP_OKAY)) != SCIP_OKAY", taking false branch.
(6) Event if_end: End of if statement.
1040 	      SCIP_CALL( SCIPallocCleanBufferArray(scip, &temp1vals, dualvectorlen) );
(7) Event cond_false: Condition "(temp2vals = BMSallocBufferMemoryArray_call(SCIPcleanbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 8UL /* sizeof (*temp2vals) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1041)) == NULL", taking false branch.
(8) Event cond_false: Condition "(_restat_ = (((temp2vals = BMSallocBufferMemoryArray_call(SCIPcleanbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 8UL /* sizeof (*temp2vals) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1041)) == NULL) ? SCIP_NOMEMORY : SCIP_OKAY)) != SCIP_OKAY", taking false branch.
(9) Event if_end: End of if statement.
1041 	      SCIP_CALL( SCIPallocCleanBufferArray(scip, &temp2vals, dualvectorlen) );
(10) Event cond_false: Condition "(temp1inds = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 4UL /* sizeof (*temp1inds) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1042)) == NULL", taking false branch.
(11) Event cond_false: Condition "(_restat_ = (((temp1inds = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 4UL /* sizeof (*temp1inds) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1042)) == NULL) ? SCIP_NOMEMORY : SCIP_OKAY)) != SCIP_OKAY", taking false branch.
(12) Event if_end: End of if statement.
1042 	      SCIP_CALL( SCIPallocBufferArray(scip, &temp1inds, dualvectorlen) );
(13) Event cond_false: Condition "(temp2inds = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 4UL /* sizeof (*temp2inds) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1043)) == NULL", taking false branch.
(14) Event cond_false: Condition "(_restat_ = (((temp2inds = BMSallocBufferMemoryArray_call(SCIPbuffer(scip), (size_t)(ptrdiff_t)dualvectorlen, 4UL /* sizeof (*temp2inds) */, "/sapmnt/home1/d029903/my_work/SCIPSoPlex_coverity/scip-soplex-9.0.0_rc1/scip/src/scip/sepa_lagromory.c", 1043)) == NULL) ? SCIP_NOMEMORY : SCIP_OKAY)) != SCIP_OKAY", taking false branch.
(15) Event if_end: End of if statement.
1043 	      SCIP_CALL( SCIPallocBufferArray(scip, &temp2inds, dualvectorlen) );
(16) Event cond_true: Condition "i < dualvectorlen", taking true branch.
(18) Event loop_begin: Jumped back to beginning of loop.
(19) Event cond_false: Condition "i < dualvectorlen", taking false branch.
1044 	      for( int i = 0; i < dualvectorlen; i++ )
1045 	      {
1046 	         temp1inds[i] = -1;
1047 	         temp2inds[i] = -1;
(17) Event loop: Jumping back to the beginning of the loop.
(20) Event loop_end: Reached end of loop.
1048 	      }
1049 	      temp2len = 0;
1050 	
1051 	      temp1vals[0] = REALABS(dualvector[0]);
1052 	      temp1inds[0] = 0;
(21) Event assignment: Assigning: "temp1len" = "1".
Also see events: [assignment][incr][divide_by_zero]
1053 	      temp1len = 1;
1054 	      pivotparam = REALABS(dualvector[0]) - radius;
1055 	
(22) Event cond_false: Condition "i < dualvectorlen", taking false branch.
1056 	      for( int i = 1; i < dualvectorlen; i++ )
1057 	      {
1058 	         if( SCIPisGT(scip, REALABS(dualvector[i]), pivotparam) )
1059 	         {
1060 	            pivotparam += ((REALABS(dualvector[i]) - pivotparam) / (temp1len + 1));
1061 	            if( SCIPisGT(scip, pivotparam, REALABS(dualvector[i]) - radius) )
1062 	            {
1063 	               temp1vals[temp1len] = REALABS(dualvector[i]);
1064 	               temp1inds[temp1len] = i;
1065 	               temp1len++;
1066 	            }
1067 	            else
1068 	            {
1069 	               for( int j = 0; j < temp1len; j++ )
1070 	               {
1071 	                  temp2vals[temp2len + j] = temp1vals[j];
1072 	                  temp2inds[temp2len + j] = temp1inds[j];
1073 	               }
1074 	               temp2len += temp1len;
1075 	               temp1vals[0] = REALABS(dualvector[i]);
1076 	               temp1inds[0] = i;
1077 	               temp1len = 1;
1078 	               pivotparam = REALABS(dualvector[i]) - radius;
1079 	            }
1080 	         }
(23) Event loop_end: Reached end of loop.
1081 	      }
1082 	
(24) Event cond_false: Condition "i < temp2len", taking false branch.
1083 	      for( int i = 0; i < temp2len; i++ )
1084 	      {
1085 	         if( SCIPisGT(scip, temp2vals[i], pivotparam) )
1086 	         {
1087 	            temp1vals[temp1len] = temp2vals[i];
1088 	            temp1inds[temp1len] = temp2inds[i];
1089 	            temp1len++;
1090 	            pivotparam += ((temp2vals[i] - pivotparam) / temp1len);
1091 	         }
(25) Event loop_end: Reached end of loop.
1092 	      }
1093 	
1094 	      temp1changed = TRUE;
(26) Event assignment: Assigning: "ntemp1removed" = "0".
Also see events: [assignment][incr][divide_by_zero]
1095 	      ntemp1removed = 0;
(27) Event cond_true: Condition "temp1changed", taking true branch.
1096 	      while( temp1changed )
1097 	      {
1098 	         temp1changed = FALSE;
1099 	
(28) Event cond_true: Condition "i < temp1len", taking true branch.
1100 	         for( int i = 0; i < temp1len; i++ )
1101 	         {
(29) Event cond_true: Condition "temp1inds[i] >= 0", taking true branch.
(30) Event cond_true: Condition "temp1vals[i] - pivotparam <= scip->set->num_epsilon", taking true branch.
1102 	            if( (temp1inds[i] >= 0) && SCIPisLE(scip, temp1vals[i], pivotparam) )
1103 	            {
1104 	               temp1inds[i] = -1;
1105 	               temp1changed = TRUE;
(31) Event incr: Incrementing "ntemp1removed". The value of "ntemp1removed" is now 1.
Also see events: [assignment][assignment][divide_by_zero]
1106 	               ntemp1removed++;
1107 	               assert(temp1len - ntemp1removed > 0);
(32) Event divide_by_zero: In expression "(pivotparam - temp1vals[i]) / (temp1len - ntemp1removed)", division by expression "temp1len - ntemp1removed" which may be zero has undefined behavior.
Also see events: [assignment][assignment][incr]
1108 	               pivotparam += ((pivotparam - temp1vals[i]) / (temp1len - ntemp1removed));
1109 	            }
1110 	         }
1111 	      }
1112 	
1113 	      for( int i = 0; i < dualvectorlen; i++ )
1114 	      {
1115 	         term = REALABS(dualvector[i]);
1116 	         val = MAX(term - pivotparam, 0.0);
1117 	
1118 	         if( SCIPisPositive(scip, dualvector[i]) )
1119 	         {
1120 	            dualvector[i] = val;
1121 	         }
1122 	         else if( SCIPisNegative(scip, dualvector[i]) )
1123 	         {
1124 	            dualvector[i] = -val;
1125 	         }
1126 	      }
1127 	
1128 	      /* free memory */
1129 	      for( int i = 0; i < dualvectorlen; i++ )
1130 	      {
1131 	         temp2vals[i] = 0.0;
1132 	         temp1vals[i] = 0.0;
1133 	      }
1134 	      SCIPfreeBufferArray(scip, &temp2inds);
1135 	      SCIPfreeBufferArray(scip, &temp1inds);
1136 	      SCIPfreeCleanBufferArray(scip, &temp2vals);
1137 	      SCIPfreeCleanBufferArray(scip, &temp1vals);
1138 	   }
1139 	
1140 	   return SCIP_OKAY;
1141 	}
1142 	
1143 	/** projection of Lagrangian multipliers onto L2-norm ball */
1144 	static
1145 	void l2BallProjection(
1146 	   SCIP*                 scip,               /**< SCIP data structure */
1147 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers to be projected onto L2-norm vall */
1148 	   int                   dualvectorlen,      /**< length of the Lagrangian multipliers vector */
1149 	   SCIP_Real             radius              /**< radius of the L2-norm ball */
1150 	   )
1151 	{
1152 	   SCIP_Real l2norm;
1153 	   SCIP_Real factor;
1154 	
1155 	   assert(!SCIPisNegative(scip, radius));
1156 	
1157 	   l2norm = 0.0;
1158 	   /* calculate the L2-norm of the Lagrangian multipliers */
1159 	   for( int i = 0; i < dualvectorlen; i++ )
1160 	   {
1161 	      l2norm += SQR(dualvector[i]);
1162 	   }
1163 	   l2norm = SQRT(l2norm);
1164 	   factor = radius/(1.0 + l2norm);
1165 	
1166 	   /* if the vector of Lagrangian multipliers is outside the L2-norm ball, then do the projection */
1167 	   if( SCIPisGT(scip, l2norm, radius) && SCIPisLT(scip, factor, 1.0) )
1168 	   {
1169 	      for( int i = 0; i < dualvectorlen; i++ )
1170 	      {
1171 	         dualvector[i] *= factor;
1172 	      }
1173 	   }
1174 	}
1175 	
1176 	/** projection of Lagrangian multipliers onto L_infinity-norm ball */
1177 	static
1178 	void linfBallProjection(
1179 	   SCIP*                 scip,               /**< SCIP data structure */
1180 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers to be projected onto L_infinity-norm vall */
1181 	   int                   dualvectorlen,      /**< length of the Lagrangian multipliers vector */
1182 	   SCIP_Real             radius              /**< radius of the L_infinity-norm ball */
1183 	   )
1184 	{
1185 	   assert(!SCIPisNegative(scip, radius));
1186 	
1187 	   /* if the vector of Lagrangian multipliers is outside the L_infinity-norm ball, then do the projection */
1188 	   for( int i = 0; i < dualvectorlen; i++ )
1189 	   {
1190 	      if( SCIPisLT(scip, dualvector[i], -radius) )
1191 	      {
1192 	         dualvector[i] = -radius;
1193 	      }
1194 	      else if( SCIPisGT(scip, dualvector[i], radius) )
1195 	      {
1196 	         dualvector[i] = radius;
1197 	      }
1198 	   }
1199 	}
1200 	
1201 	/** weighted Lagrangian multipliers based on a given vector as stability center */
1202 	/* @todo calculate weight outside this function and pass it (so that this function becomes generic and independent of
1203 	 * the terminology related to best Lagrangian multipliers)
1204 	 */
1205 	static
1206 	SCIP_RETCODE weightedDualVector(
1207 	   SCIP*                 scip,               /**< SCIP data structure */
1208 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1209 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers */
1210 	   int                   dualvectorlen,      /**< length of the Lagrangian multipliers vector */
1211 	   SCIP_Real*            stabilitycenter,    /**< stability center (i.e., core vector of Lagrangian multipliers) */
1212 	   int                   stabilitycenterlen, /**< length of the stability center */
1213 	   int                   nbestdualupdates,   /**< number of best Lagrangian values found so far */
1214 	   int                   totaliternum        /**< total number of iterations of the relax-and-cut algorithm performed so far */
1215 	   )
1216 	{
1217 	   SCIP_Real constant;
1218 	   SCIP_Real weight;
1219 	   SCIP_Real alpha;
1220 	
1221 	   constant = MAX(2.0, sepadata->constant);
1222 	   /* weight factor from the literature on Dantzig-Wolfe decomposition stabilization schemes */
1223 	   weight = MIN(constant, (totaliternum + 1 + nbestdualupdates) / 2.0);
1224 	   alpha = 1.0 / weight;
1225 	
1226 	   assert(dualvectorlen >= stabilitycenterlen);
1227 	
1228 	   /* weighted Lagrangian multipliers */
1229 	   for( int i = 0; i < stabilitycenterlen; i++ )
1230 	   {
1231 	      dualvector[i] = alpha * dualvector[i] + (1 - alpha) * stabilitycenter[i];
1232 	   }
1233 	   for( int i = stabilitycenterlen; i < dualvectorlen; i++ )
1234 	   {
1235 	      dualvector[i] = alpha * dualvector[i];
1236 	   }
1237 	
1238 	   return SCIP_OKAY;
1239 	}
1240 	
1241 	/** stabilize Lagrangian multipliers */
1242 	static
1243 	SCIP_RETCODE stabilizeDualVector(
1244 	   SCIP*                 scip,               /**< SCIP data structure */
1245 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1246 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers */
1247 	   int                   dualvectorlen,      /**< length of the Lagrangian multipliers vector */
1248 	   SCIP_Real*            bestdualvector,     /**< best Lagrangian multipliers found so far */
1249 	   int                   bestdualvectorlen,  /**< length of the best Lagrangian multipliers vector */
1250 	   int                   nbestdualupdates,   /**< number of best Lagrangian values found so far */
1251 	   int                   subgradientiternum, /**< iteration number of the subgradient algorithm */
1252 	   int                   totaliternum,       /**< total number of iterations of the relax-and-cut algorithm performed so far */
1253 	   SCIP_Real             maxviolscore,       /**< weighted average of maximum value of generated cut violations and maximum value of
1254 	                                               complementarity slackness violations, in the current iteration */
1255 	   SCIP_Real             maxviolscoreold,    /**< weighted average of maximum value of generated cut violations and maximum value of
1256 	                                               complementarity slackness violations, in the previous iteration */
1257 	   SCIP_Real             nviolscore,         /**< weighted average of number of generated cut violations and number of complementarity
1258 	                                               slackness violations, in the current iteration */
1259 	   SCIP_Real             nviolscoreold,      /**< weighted average of number of generated cut violations and number of complementarity
1260 	                                               slackness violations, in the previous iteration */
1261 	   int                   nlpiters,           /**< number of LP iterations taken for solving the Lagrangian dual with fixed multipliers in
1262 	                                               current iteration */
1263 	   SCIP_Real*            ballradius          /**< norm ball radius */
1264 	   )
1265 	{
1266 	   if( sepadata->projectiontype > 0 )
1267 	   {
1268 	      if( subgradientiternum >= 1 )
1269 	      {
1270 	         /* update the ball radius */
1271 	         updateBallRadius(scip, sepadata, maxviolscore, maxviolscoreold, nviolscore, nviolscoreold, nlpiters,
1272 	               ballradius);
1273 	      }
1274 	
1275 	      if( sepadata->projectiontype == 1 )
1276 	      {
1277 	         /* projection of Lagrangian multipliers onto L1-norm ball */
1278 	         SCIP_CALL( l1BallProjection(scip, dualvector, dualvectorlen, *ballradius) );
1279 	      }
1280 	      else if( sepadata->projectiontype == 2 )
1281 	      {
1282 	         /* projection of Lagrangian multipliers onto L2-norm ball */
1283 	         l2BallProjection(scip, dualvector, dualvectorlen, *ballradius);
1284 	      }
1285 	      else if( sepadata->projectiontype == 3 )
1286 	      {
1287 	         /* projection of Lagrangian multipliers onto L_inf-norm ball */
1288 	         linfBallProjection(scip, dualvector, dualvectorlen, *ballradius);
1289 	      }
1290 	   }
1291 	
1292 	   if( sepadata->stabilitycentertype == 1 )
1293 	   {
1294 	      /* weighted Lagrangian multipliers based on best Langrangian multipliers as stability center */
1295 	      SCIP_CALL( weightedDualVector(scip, sepadata, dualvector, dualvectorlen, bestdualvector,
1296 	               bestdualvectorlen, nbestdualupdates, totaliternum) );
1297 	   }
1298 	
1299 	   return SCIP_OKAY;
1300 	}
1301 	
1302 	/** update Lagrangian multipliers */
1303 	static
1304 	SCIP_RETCODE updateDualVector(
1305 	   SCIP*                 scip,               /**< SCIP data structure */
1306 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1307 	   SCIP_Real*            dualvector1,        /**< Lagrangian multipliers vector to be updated */
1308 	   SCIP_Real*            dualvector2,        /**< Lagrangian multipliers vector used for backtracking */
1309 	   int                   dualvector2len,     /**< length of the Lagrangian multipliers vector used for backtracking */
1310 	   int                   ndualvector2updates,/**< number of best Lagrangian values found so far */
1311 	   int                   subgradientiternum, /**< iteration number of the subgradient algorithm */
1312 	   int                   totaliternum,       /**< total number of iterations of the relax-and-cut algorithm performed so far */
1313 	   SCIP_Real             steplength,         /**< step length used for updating Lagrangian multipliers */
1314 	   SCIP_Real*            subgradient,        /**< subgradient vector */
1315 	   int                   ncuts,              /**< number of generated cuts so far */
1316 	   SCIP_Bool             backtrack,          /**< whether the Lagrangian multipliers need to be backtracked */
1317 	   SCIP_Real             maxviolscore,       /**< weighted average of maximum value of generated cut violations and maximum value of
1318 	                                               complementarity slackness violations, in the current iteration */
1319 	   SCIP_Real             maxviolscoreold,    /**< weighted average of maximum value of generated cut violations and maximum value of
1320 	                                               complementarity slackness violations, in the previous iteration */
1321 	   SCIP_Real             nviolscore,         /**< weighted average of number of generated cut violations and number of complementarity
1322 	                                               slackness violations, in the current iteration */
1323 	   SCIP_Real             nviolscoreold,      /**< weighted average of number of generated cut violations and number of complementarity
1324 	                                               slackness violations, in the previous iteration */
1325 	   int                   nlpiters,           /**< number of LP iterations taken for solving the Lagrangian dual with fixed multipliers in
1326 	                                               current iteration */
1327 	   SCIP_Bool*            dualvecsdiffer,     /**< whether the updated Lagrangian multipliers differ from the old one */
1328 	   SCIP_Real*            ballradius          /**< norm ball radius */
1329 	   )
1330 	{
1331 	   SCIP_Real* dualvector1copy;
1332 	
1333 	   assert(dualvector2len <= ncuts);
1334 	   assert(dualvecsdiffer != NULL);
1335 	
1336 	   *dualvecsdiffer = FALSE;
1337 	   /* @todo do allocation and free operations outside only once instead of every time this function is called? */
1338 	   /* copy of the Lagrangian multipliers to be used to check if the updated vector is different than this */
1339 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &dualvector1copy, ncuts) );
1340 	   for( int i = 0; i < ncuts; i++ )
1341 	   {
1342 	      dualvector1copy[i] = dualvector1[i];
1343 	   }
1344 	
1345 	   /* if backtracking was not identified at the time of the mu parameter update, then update the Lagrangian multipliers
1346 	    * based on the given subgradient vector
1347 	    */
1348 	   if( !backtrack )
1349 	   {
1350 	      assert((subgradient != NULL) || (ncuts == 0));
1351 	      assert(subgradientiternum >= 0);
1352 	
1353 	      /* update Lagrangian multipliers */
1354 	      for( int i = 0; i < ncuts; i++ )
1355 	      {
1356 	         dualvector1[i] += steplength * subgradient[i];
1357 	      }
1358 	
1359 	      /* projection onto non-negative orthant */
1360 	      for( int i = 0; i < ncuts; i++ )
1361 	      {
1362 	         dualvector1[i] = MAX(dualvector1[i], 0.0);
1363 	      }
1364 	
1365 	      /* stabilization of Lagrangian multipliers */
1366 	      SCIP_CALL( stabilizeDualVector(scip, sepadata, dualvector1, ncuts, dualvector2, dualvector2len,
1367 	               ndualvector2updates, subgradientiternum, totaliternum, maxviolscore, maxviolscoreold, nviolscore,
1368 	               nviolscoreold, nlpiters, ballradius) );
1369 	
1370 	      /* projection onto non-negative orthant again in case stabilization changed some components negative*/
1371 	      for( int i = 0; i < ncuts; i++ )
1372 	      {
1373 	         dualvector1[i] = MAX(dualvector1[i], 0.0);
1374 	      }
1375 	   }
1376 	   /* if backtracking was identified at the time of the mu parameter update, then backtrack the Lagrangian multipliers
1377 	    * based on the given backtracking multipliers
1378 	    */
1379 	   else
1380 	   {
1381 	      for( int i = 0; i < dualvector2len; i++ )
1382 	      {
1383 	         dualvector1[i] = dualvector2[i];
1384 	      }
1385 	
1386 	      for( int i = dualvector2len; i < ncuts; i++ )
1387 	      {
1388 	         dualvector1[i] = 0.0;
1389 	      }
1390 	   }
1391 	
1392 	   /* identify if the vector of Lagrangian multipliers is indeed different after updating */
1393 	   for( int i = 0; i < ncuts; i++ )
1394 	   {
1395 	      if( !SCIPisEQ(scip, dualvector1[i], dualvector1copy[i]) )
1396 	      {
1397 	         *dualvecsdiffer = TRUE;
1398 	         break;
1399 	      }
1400 	   }
1401 	
1402 	   /* free memory */
1403 	   for( int i = 0; i < ncuts; i++ )
1404 	   {
1405 	      dualvector1copy[i] = 0.0;
1406 	   }
1407 	   SCIPfreeCleanBufferArray(scip, &dualvector1copy);
1408 	
1409 	   return SCIP_OKAY;
1410 	}
1411 	
1412 	/** check different termination criteria */
1413 	/* @note: the criterion based on objvecsdiffer assumes deterministic solving process (i.e., we would get same LP solution
1414 	 * for "Lagrangian dual with fixed Lagrangian multipliers" when the objective vector remains the same across iterations).
1415 	 */
1416 	/* @todo nlpssolved criterion? */
1417 	static
1418 	SCIP_RETCODE checkLagrangianDualTermination(
1419 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
1420 	   int                   nnewaddedsoftcuts,        /**< number of cuts that were recently penalized and added to the Lagrangian dual's
1421 	                                                     objective function */
1422 	   int                   nyettoaddsoftcuts,        /**< number of cuts that are yet to be penalized and added to the Lagrangian dual's
1423 	                                                     objective function */
1424 	   SCIP_Bool             objvecsdiffer,            /**< whether the Lagrangian dual's objective function has changed */
1425 	   int                   ngeneratedcurrroundcuts,  /**< number of cuts generated in the current separation round */
1426 	   int                   nmaxgeneratedperroundcuts,/**< maximal number of cuts allowed to generate per separation round */
1427 	   int                   ncurrroundlpiters,        /**< number of separating LP iterations in the current separation round */
1428 	   int                   depth,                    /**< depth of the current node */
1429 	   SCIP_Bool*            terminate                 /**< whether to terminate the subgradient algorithm loop */
1430 	   )
1431 	{
1432 	   *terminate = FALSE;
1433 	
1434 	   /* check if no new cuts were added to the Lagrangian dual, no cuts are remaining to be added, and the objective
1435 	    * function of the Lagrangian dual with fixed multipliers had not changed from the previous iteration
1436 	    */
1437 	   if( (nnewaddedsoftcuts == 0) && (nyettoaddsoftcuts == 0) && !objvecsdiffer )
1438 	      *terminate = TRUE;
1439 	
1440 	   /* check if allowed number of cuts in this separation round have already been generated */
1441 	   if( ngeneratedcurrroundcuts >= nmaxgeneratedperroundcuts )
1442 	      *terminate = TRUE;
1443 	
1444 	   /* check if allowed number of cuts in the tree have already been generated */
1445 	   if( sepadata->ntotalcuts >= sepadata->nmaxtotalcuts )
1446 	      *terminate = TRUE;
1447 	
1448 	   /* check if allowed number of simplex iterations in this separation round have already been used up */
1449 	   if( (sepadata->nmaxperroundlpiters >= 0) && (ncurrroundlpiters >= sepadata->nmaxperroundlpiters) )
1450 	      *terminate = TRUE;
1451 	
1452 	   /* check if allowed number of simplex iterations in the root node have already been used up */
1453 	   if( (depth == 0) && (sepadata->nmaxrootlpiters >= 0) && (sepadata->nrootlpiters >= sepadata->nmaxrootlpiters) )
1454 	      *terminate = TRUE;
1455 	
1456 	   /* check if allowed number of simplex iterations in the tree have already been used up */
1457 	   if( (sepadata->nmaxtotallpiters >= 0) && (sepadata->ntotallpiters >= sepadata->nmaxtotallpiters) )
1458 	      *terminate = TRUE;
1459 	
1460 	   return SCIP_OKAY;
1461 	}
1462 	
1463 	/** solve the LP corresponding to the Lagrangian dual with fixed Lagrangian multipliers */
1464 	static
1465 	SCIP_RETCODE solveLagromoryLP(
1466 	   SCIP*                 scip,               /**< SCIP data structure */
1467 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1468 	   int                   depth,              /**< depth of the current node in the tree */
1469 	   SCIP_Real             origobjoffset,      /**< objective offset in the current node's relaxation */
1470 	   SCIP_Bool*            solfound,           /**< whether an LP optimal solution has been found */
1471 	   SCIP_SOL*             sol,                /**< data structure to store LP optimal solution, if found */
1472 	   SCIP_Real*            solvals,            /**< values of the LP optimal solution, if found */
1473 	   SCIP_Real*            objval,             /**< optimal objective value of the LP optimal solution, if found */
1474 	   int*                  ncurrroundlpiters   /**< number of LP iterations taken for solving Lagrangian dual problems with fixed multipliers
1475 	                                               in the current separator round */
1476 	   )
1477 	{
1478 	   SCIP_Real timelimit;
1479 	   SCIP_COL** cols;
1480 	   SCIP_COL* col;
1481 	   SCIP_VAR* var;
1482 	   SCIP_LPI* lpi;
1483 	   SCIP_Bool lperror;
1484 	   SCIP_Bool cutoff;
1485 	   SCIP_LPSOLSTAT stat;
1486 	   SCIP_Longint ntotallpiters;
1487 	   SCIP_Longint nlpiters;
1488 	   int ncols;
1489 	   int iterlimit;
1490 	
1491 	   assert(solfound != NULL);
1492 	   assert(sol != NULL);
1493 	   assert(solvals != NULL);
1494 	   assert(ncurrroundlpiters != NULL);
1495 	   assert(objval != NULL);
1496 	
1497 	   *solfound = FALSE;
1498 	   lperror = FALSE;
1499 	   cutoff = FALSE;
1500 	   iterlimit = -1;
1501 	   lpi = sepadata->lpiwithsoftcuts;
1502 	
1503 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
1504 	
1505 	   /* set time limit */
1506 	   SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1507 	   if( !SCIPisInfinity(scip, timelimit) )
1508 	   {
1509 	      timelimit -= SCIPgetSolvingTime(scip);
1510 	      if( timelimit <= 0.0 )
1511 	      {
1512 	         SCIPdebugMsg(scip, "skip Lagromory cut generation since no time left\n");
1513 	         goto TERMINATE;
1514 	      }
1515 	      /* @note: the following direct LPI call is being used because of the lack of an equivalent function call in
1516 	       * scip_lp.c (lpSetRealpar exists in lp.c though)
1517 	       */
1518 	      SCIP_CALL( SCIPlpiSetRealpar(lpi, SCIP_LPPAR_LPTILIM, timelimit) );
1519 	   }
1520 	
1521 	   /* find iteration limit */
1522 	   if( (depth == 0) &&
1523 	         (sepadata->perrootlpiterfactor >= 0.0 && !SCIPisInfinity(scip, sepadata->perrootlpiterfactor)) )
1524 	   {
1525 	      iterlimit = (int)(sepadata->perrootlpiterfactor * SCIPgetNRootFirstLPIterations(scip));
1526 	   }
1527 	   else if( (depth > 0) &&
1528 	         (sepadata->perlpiterfactor >= 0.0 && !SCIPisInfinity(scip, sepadata->perlpiterfactor)) )
1529 	   {
1530 	      iterlimit = (int)(sepadata->perlpiterfactor * SCIPgetNNodeInitLPIterations(scip));
1531 	   }
1532 	   if( sepadata->nmaxperroundlpiters >= 0 )
1533 	   {
1534 	      if( sepadata->nmaxperroundlpiters - *ncurrroundlpiters >= 0 )
1535 	      {
1536 	         if( iterlimit >= 0 )
1537 	         {
1538 	            iterlimit = MIN(iterlimit, sepadata->nmaxperroundlpiters - *ncurrroundlpiters);
1539 	         }
1540 	         else
1541 	         {
1542 	            iterlimit = sepadata->nmaxperroundlpiters - *ncurrroundlpiters;
1543 	         }
1544 	      }
1545 	      else
1546 	      {
1547 	         iterlimit = 0;
1548 	      }
1549 	   }
1550 	   /* @todo impose a finite iteration limit only when the dualvector changes from zero to non-zero for the first time because
1551 	    * many simplex pivots are performed in this case even with warm starting (compared to the case when the
1552 	    * dualvector changes from non-zero to non-zero).
1553 	    */
1554 	
1555 	   /* solve the LP with an iteration limit and get number of simplex iterations taken */
1556 	   ntotallpiters = SCIPgetNLPIterations(scip);
1557 	
1558 	   SCIP_CALL( SCIPsolveDiveLP(scip, iterlimit, &lperror, &cutoff) );
1559 	
1560 	   nlpiters = SCIPgetNLPIterations(scip) - ntotallpiters;
1561 	
1562 	   /* get the solution and objective value if optimal */
1563 	   stat = SCIPgetLPSolstat(scip);
1564 	   /* @todo is there any way to accept terminations due to iterlimit and timelimit as well? It is not possible
1565 	    * currently because primal sol is not saved in these cases.
1566 	    */
1567 	   /* @note: ideally, only primal feasibility is sufficient. But, there is no such option with SCIPgetLPSolstat. */
1568 	   if( stat == SCIP_LPSOLSTAT_OPTIMAL )
1569 	   {
1570 	      if( SCIPisLPSolBasic(scip) )
1571 	      {
1572 	         *solfound = TRUE;
1573 	
1574 	         /* update sol */
1575 	         for( int i = 0; i < ncols; ++i )
1576 	         {
1577 	            col = cols[i];
1578 	            assert(col != NULL);
1579 	
1580 	            var = SCIPcolGetVar(col);
1581 	
1582 	            solvals[i] = SCIPcolGetPrimsol(col);
1583 	            SCIP_CALL( SCIPsetSolVal(scip, sol, var, solvals[i]) );
1584 	         }
1585 	
1586 	         *objval = SCIPgetLPObjval(scip);
1587 	         *objval = *objval + origobjoffset;
1588 	      }
1589 	   }
1590 	
1591 	   /* update some statistics */
1592 	   if( depth == 0 )
1593 	   {
1594 	      sepadata->nrootlpiters += (int)nlpiters;
1595 	   }
1596 	   sepadata->ntotallpiters += (int)nlpiters;
1597 	   *ncurrroundlpiters += (int)nlpiters;
1598 	
1599 	TERMINATE:
1600 	   return SCIP_OKAY;
1601 	}
1602 	
1603 	/** solve the LP corresponding to the node relaxation upon adding all the generated cuts */
1604 	static
1605 	SCIP_RETCODE solveLPWithHardCuts(
1606 	   SCIP*                 scip,               /**< SCIP data structure */
1607 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1608 	   SCIP_Bool*            solfound,           /**< whether an LP optimal solution has been found */
1609 	   SCIP_SOL*             sol,                /**< data structure to store LP optimal solution, if found */
1610 	   SCIP_Real*            solvals             /**< values of the LP optimal solution, if found */
1611 	   )
1612 	{
1613 	   SCIP_Real timelimit;
1614 	   SCIP_COL** cols;
1615 	   SCIP_COL* col;
1616 	   SCIP_VAR* var;
1617 	   int ncols;
1618 	
1619 	   assert(solfound != NULL);
1620 	   assert(sol != NULL);
1621 	   assert(solvals != NULL);
1622 	
1623 	   *solfound = FALSE;
1624 	
1625 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
1626 	
1627 	   /* set time limit */
1628 	   SCIP_CALL( SCIPgetRealParam(scip, "limits/time", &timelimit) );
1629 	   if( !SCIPisInfinity(scip, timelimit) )
1630 	   {
1631 	      timelimit -= SCIPgetSolvingTime(scip);
1632 	      if( timelimit <= 0.0 )
1633 	      {
1634 	         SCIPdebugMsg(scip, "skip Lagromory cut generation since no time left\n");
1635 	         goto TERMINATE;
1636 	      }
1637 	      SCIP_CALL( SCIPlpiSetRealpar(sepadata->lpiwithhardcuts, SCIP_LPPAR_LPTILIM, timelimit) );
1638 	   }
1639 	
1640 	   /* solve the LP */
1641 	   SCIPlpiSolvePrimal(sepadata->lpiwithhardcuts);  /*lint !e534*/
1642 	
1643 	   /* get the solution if primal feasible */
1644 	   if( SCIPlpiIsPrimalFeasible(sepadata->lpiwithhardcuts) )
1645 	   {
1646 	      *solfound = TRUE;
1647 	      SCIP_CALL( SCIPlpiGetSol(sepadata->lpiwithhardcuts, NULL, solvals, NULL, NULL, NULL) );
1648 	
1649 	      /* update sol */
1650 	      for( int i = 0; i < ncols; ++i )
1651 	      {
1652 	         col = cols[i];
1653 	         assert(col != NULL);
1654 	
1655 	         var = SCIPcolGetVar(col);
1656 	
1657 	         SCIP_CALL( SCIPsetSolVal(scip, sol, var, solvals[i]) );
1658 	      }
1659 	   }
1660 	
1661 	TERMINATE:
1662 	   return SCIP_OKAY;
1663 	}
1664 	
1665 	/** construct a cut based on the input cut coefficients, sides, etc */
1666 	static
1667 	SCIP_RETCODE constructCutRow(
1668 	   SCIP*                 scip,                     /**< SCIP data structure */
1669 	   SCIP_SEPA*            sepa,                     /**< pointer to the separator */
1670 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
1671 	   int                   mainiternum,              /**< iteration number of the outer loop of the relax-and-cut algorithm */
1672 	   int                   subgradientiternum,       /**< iteration number of the subgradient algorithm */
1673 	   int                   cutnnz,                   /**< number of nonzeros in cut */
1674 	   int*                  cutinds,                  /**< column indices in cut */
1675 	   SCIP_Real*            cutcoefs,                 /**< cut cofficients */
1676 	   SCIP_Real             cutefficacy,              /**< cut efficacy */
1677 	   SCIP_Real             cutrhs,                   /**< RHS of cut */
1678 	   SCIP_Bool             cutislocal,               /**< whether cut is local */
1679 	   int                   cutrank,                  /**< rank of cut */
1680 	   SCIP_ROW**            generatedcuts,            /**< array of generated cuts */
1681 	   SCIP_Real*            generatedcutefficacies,   /**< array of generated cut efficacies w.r.t. the respective LP bases used for cut
1682 	                                                     generations */
1683 	   int                   ngeneratedcurrroundcuts,  /**< number of cuts generated until the previous basis in the current separation round */
1684 	   int*                  ngeneratednewcuts,        /**< number of new cuts generated using the current basis */
1685 	   SCIP_Bool*            cutoff                    /**< should the current node be cutoff? */
1686 	   )
1687 	{
1688 	   SCIP_COL** cols;
1689 	   SCIP_Real minactivity;
1690 	   SCIP_Real maxactivity;
1691 	   SCIP_Real lhs;
1692 	   SCIP_Real rhs;
1693 	
1694 	   assert(scip != NULL);
1695 	   assert(mainiternum >= 0);
1696 	   assert(ngeneratednewcuts != NULL);
1697 	   assert(cutoff != NULL);
1698 	
1699 	   *cutoff = FALSE;
1700 	
1701 	   if( cutnnz == 0 && SCIPisFeasNegative(scip, cutrhs) ) /*lint !e644*/
1702 	   {
1703 	      SCIPdebugMsg(scip, " -> Lagromory cut detected node infeasibility with cut 0 <= %g.\n", cutrhs);
1704 	      *cutoff = TRUE;
1705 	      return SCIP_OKAY;
1706 	   }
1707 	
1708 	   /* only take efficient cuts */
1709 	   if( SCIPisEfficacious(scip, cutefficacy) )
1710 	   {
1711 	      SCIP_ROW* cut;
1712 	      SCIP_VAR* var;
1713 	      char cutname[SCIP_MAXSTRLEN];
1714 	      int v;
1715 	
1716 	      /* construct cut name */
1717 	      if( subgradientiternum >= 0 )
1718 	      {
1719 	         (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "%s_%" SCIP_LONGINT_FORMAT "_%d" "_%d" "_%d", SCIPsepaGetName(sepa),
1720 	               sepadata->ncalls, mainiternum, subgradientiternum, ngeneratedcurrroundcuts + *ngeneratednewcuts);
1721 	      }
1722 	      else
1723 	      {
1724 	         (void) SCIPsnprintf(cutname, SCIP_MAXSTRLEN, "%s_%" SCIP_LONGINT_FORMAT "_%d" "_%d", SCIPsepaGetName(sepa),
1725 	               sepadata->ncalls, mainiternum, ngeneratedcurrroundcuts + *ngeneratednewcuts);
1726 	      }
1727 	
1728 	      /* create empty cut */
1729 	      SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &cut, sepa, cutname, -SCIPinfinity(scip), cutrhs, cutislocal, FALSE,
1730 	               sepadata->dynamiccuts) );
1731 	
1732 	      /* set cut rank */
1733 	      SCIProwChgRank(cut, cutrank); /*lint !e644*/
1734 	
1735 	      /* cache the row extension and only flush them if the cut gets added */
1736 	      SCIP_CALL( SCIPcacheRowExtensions(scip, cut) );
1737 	
1738 	      cols = SCIPgetLPCols(scip);
1739 	
1740 	      /* collect all non-zero coefficients */
1741 	      for( v = 0; v < cutnnz; ++v )
1742 	      {
1743 	         var = SCIPcolGetVar(cols[cutinds[v]]);
1744 	         SCIP_CALL( SCIPaddVarToRow(scip, cut, var, cutcoefs[v]) );
1745 	      }
1746 	
1747 	      /* flush all changes before adding the cut */
1748 	      SCIP_CALL( SCIPflushRowExtensions(scip, cut) );
1749 	
1750 	      if( SCIProwGetNNonz(cut) == 0 )
1751 	      {
1752 	         assert( SCIPisFeasNegative(scip, cutrhs) );
1753 	         SCIPdebugMsg(scip, " -> Lagromory cut detected node infeasibility with cut 0 <= %g.\n", cutrhs);
1754 	         *cutoff = TRUE;
1755 	         return SCIP_OKAY;
1756 	      }
1757 	      else
1758 	      {
1759 	         /* gathering lhs and rhs again in case the separator is extended later to add other cuts/constraints that may
1760 	          * have non-inf lhs or inf rhs */
1761 	         lhs = SCIProwGetLhs(cut);
1762 	         rhs = SCIProwGetRhs(cut);
1763 	         assert(SCIPisInfinity(scip, -lhs));
1764 	         assert(!SCIPisInfinity(scip, rhs));
1765 	
1766 	         SCIPdebugMsg(scip, " -> %s cut <%s>: rhs=%f, eff=%f\n", "lagromory", cutname, cutrhs, cutefficacy);
1767 	
1768 	         /* check if the cut leads to infeasibility (i.e., *cutoff = TRUE) */
1769 	         {
1770 	            /* modifiable cuts cannot be declared infeasible, since we don't know all coefficients */
1771 	            if( SCIProwIsModifiable(cut) )
1772 	               *cutoff = FALSE;
1773 	
1774 	            /* check for activity infeasibility */
1775 	            minactivity = SCIPgetRowMinActivity(scip, cut);
1776 	            maxactivity = SCIPgetRowMaxActivity(scip, cut);
1777 	
1778 	            if( (!SCIPisInfinity(scip,  rhs) && SCIPisFeasPositive(scip, minactivity - rhs)) ||
1779 	                  (!SCIPisInfinity(scip, -lhs) && SCIPisFeasNegative(scip, maxactivity - lhs)) )
1780 	            {
1781 	               SCIPdebugMsg(scip, "cut <%s> is infeasible (sides=[%g,%g], act=[%g,%g])\n",
1782 	                     SCIProwGetName(cut), lhs, rhs, minactivity, maxactivity);
1783 	               *cutoff = TRUE;
1784 	            }
1785 	         }
1786 	
1787 	         /* store the newly generated cut in an array and update some statistics */
1788 	         generatedcuts[ngeneratedcurrroundcuts + *ngeneratednewcuts] = cut;
1789 	         generatedcutefficacies[ngeneratedcurrroundcuts + *ngeneratednewcuts] = cutefficacy;
1790 	         ++(*ngeneratednewcuts);
1791 	      }
1792 	   }
1793 	
1794 	   return SCIP_OKAY;
1795 	}
1796 	
1797 	/** aggregated generated cuts based on the best Lagrangian multipliers */
1798 	static
1799 	SCIP_RETCODE aggregateGeneratedCuts(
1800 	   SCIP*                 scip,               /**< SCIP data structure */
1801 	   SCIP_SEPA*            sepa,               /**< pointer to the separator */
1802 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
1803 	   SCIP_ROW**            generatedcuts,      /**< cuts generated in the current separation round */
1804 	   SCIP_Real*            bestdualvector,     /**< best Lagrangian multipliers vector */
1805 	   int                   bestdualvectorlen,  /**< length of the best Lagrangian multipliers vector */
1806 	   SCIP_ROW**            aggrcuts,           /**< aggregated cuts generated so far in the current separation round */
1807 	   int*                  naggrcuts,          /**< number of aggregated cuts generated so far in the current separation round */
1808 	   SCIP_Bool*            cutoff              /**< should the current node be cutoff? */
1809 	   )
1810 	{
1811 	   SCIP_Real* cutvals;           /**< cut cofficients */
1812 	   SCIP_COL** cutcols;
1813 	   SCIP_COL** cols;
1814 	   SCIP_VAR* var;
1815 	   SCIP_ROW* cut;
1816 	   SCIP_ROW* aggrcut;
1817 	   SCIP_Real minactivity;
1818 	   SCIP_Real maxactivity;
1819 	   SCIP_Real cutlhs;
1820 	   SCIP_Real cutrhs;
1821 	   SCIP_Real cutconst;
1822 	   SCIP_Real aggrcutlhs;
1823 	   SCIP_Real aggrcutrhs;
1824 	   SCIP_Real aggrcutconst;
1825 	   SCIP_Real* aggrcutvals;
1826 	   SCIP_Real* aggrcutcoefs;
1827 	   SCIP_Real multiplier;
1828 	   SCIP_Bool aggrcutislocal;
1829 	   SCIP_Bool aggrindicator;
1830 	   SCIP_Real* tmpcutvals;
1831 	   SCIP_Real QUAD(quadterm);
1832 	   SCIP_Real QUAD(tmpcutconst);
1833 	   SCIP_Real QUAD(tmpcutrhs);
1834 	   SCIP_Real QUAD(quadprod);
1835 	   char aggrcutname[SCIP_MAXSTRLEN];
1836 	   int cutnnz;             /**< number of nonzeros in cut */
1837 	   int aggrcutnnz;
1838 	   int* aggrcutinds;            /**< column indices in cut */
1839 	   int aggrcutrank;            /**< rank of cut */
1840 	   int cutrank;
1841 	   int ncols;
1842 	   int collppos;
1843 	   int nlocalcuts;
1844 	
1845 	   assert(scip != NULL);
1846 	   assert(naggrcuts != NULL);
1847 	   assert(cutoff != NULL);
1848 	
1849 	   *cutoff = FALSE;
1850 	   aggrcutlhs = -SCIPinfinity(scip);
1851 	   aggrcutnnz = 0;
1852 	   aggrcutrank = -1;
1853 	   nlocalcuts = 0;
1854 	   aggrindicator = FALSE;
1855 	
1856 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
1857 	
1858 	   /* allocate memory */
1859 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &tmpcutvals, QUAD_ARRAY_SIZE(ncols)) );
1860 	   QUAD_ASSIGN(tmpcutconst, 0.0);
1861 	   QUAD_ASSIGN(tmpcutrhs, 0.0);
1862 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &aggrcutvals, ncols) );
1863 	   SCIP_CALL( SCIPallocBufferArray(scip, &aggrcutcoefs, ncols) );
1864 	   SCIP_CALL( SCIPallocBufferArray(scip, &aggrcutinds, ncols) );
1865 	
1866 	   /* aggregate cuts based on the input Lagrangian multipliers */
1867 	   for( int i = 0; i < bestdualvectorlen; i++ )
1868 	   {
1869 	      multiplier = bestdualvector[i];
1870 	      if( SCIPisGE(scip, multiplier, 1e-4) )
1871 	      {
1872 	         cut = generatedcuts[i];
1873 	         cutnnz = SCIProwGetNNonz(cut);
1874 	         cutlhs = SCIProwGetLhs(cut);
1875 	         cutrhs = SCIProwGetRhs(cut);
1876 	         assert(SCIPisInfinity(scip, -cutlhs));
1877 	         assert(!SCIPisInfinity(scip, cutrhs));
1878 	         cutvals = SCIProwGetVals(cut);
1879 	         cutcols = SCIProwGetCols(cut);
1880 	         cutconst = SCIProwGetConstant(cut);
1881 	
1882 	         for( int j = 0; j < cutnnz; j++ )
1883 	         {
1884 	            collppos = SCIPcolGetLPPos(cutcols[j]);
1885 	            assert(collppos >= 0);
1886 	            assert(collppos <= ncols);
1887 	
1888 	            QUAD_ARRAY_LOAD(quadterm, tmpcutvals, collppos);
1889 	            SCIPquadprecProdDD(quadprod, multiplier, cutvals[j]);
1890 	            SCIPquadprecSumQQ(quadterm, quadterm, quadprod);
1891 	            QUAD_ARRAY_STORE(tmpcutvals, collppos, quadterm);
1892 	         }
1893 	
1894 	         SCIPquadprecProdDD(quadprod, multiplier, cutconst);
1895 	         SCIPquadprecSumQQ(tmpcutconst, tmpcutconst, quadprod);
1896 	         SCIPquadprecProdDD(quadprod, multiplier, cutrhs);
1897 	         SCIPquadprecSumQQ(tmpcutrhs, tmpcutrhs, quadprod);
1898 	
1899 	         cutrank = SCIProwGetRank(cut);
1900 	         aggrcutrank = MAX(aggrcutrank, cutrank);
1901 	         nlocalcuts += (SCIProwIsLocal(cut) ? 1 : 0);
1902 	         aggrindicator = TRUE;
1903 	      }
1904 	   }
1905 	   /* if at least one cut is local, then the aggregated cut is local */
1906 	   aggrcutislocal = (nlocalcuts > 0 ? TRUE : FALSE);
1907 	
1908 	   /* if the aggregation was successful, then create an empty row and build a cut */
1909 	   if( aggrindicator )
1910 	   {
1911 	      aggrcutconst = QUAD_TO_DBL(tmpcutconst);
1912 	      aggrcutrhs = QUAD_TO_DBL(tmpcutrhs);
1913 	
1914 	      /* build sparse representation of the aggregated cut */
1915 	      for( int i = 0; i < ncols; i++ )
1916 	      {
1917 	         QUAD_ARRAY_LOAD(quadterm, tmpcutvals, i);
1918 	         aggrcutvals[i] = QUAD_TO_DBL(quadterm);
1919 	         if( !SCIPisZero(scip, aggrcutvals[i]) )
1920 	         {
1921 	            aggrcutcoefs[aggrcutnnz] = aggrcutvals[i];
1922 	            aggrcutinds[aggrcutnnz] = i;
1923 	            aggrcutnnz++;
1924 	         }
1925 	      }
1926 	
1927 	      if( aggrcutnnz == 0 && SCIPisFeasNegative(scip, aggrcutrhs) ) /*lint !e644*/
1928 	      {
1929 	         SCIPdebugMsg(scip, " -> Lagromory cut detected node infeasibility with cut 0 <= %g.\n", aggrcutrhs);
1930 	         *cutoff = TRUE;
1931 	         goto TERMINATE;
1932 	      }
1933 	
1934 	      /* construct aggregated cut name */
1935 	      (void) SCIPsnprintf(aggrcutname, SCIP_MAXSTRLEN, "%s_%" SCIP_LONGINT_FORMAT "_aggregated", SCIPsepaGetName(sepa),
1936 	            sepadata->ncalls);
1937 	
1938 	      /* create empty cut */
1939 	      SCIP_CALL( SCIPcreateEmptyRowSepa(scip, &aggrcut, sepa, aggrcutname, aggrcutlhs - aggrcutconst, aggrcutrhs -
1940 	               aggrcutconst, aggrcutislocal, FALSE, sepadata->dynamiccuts) );
1941 	
1942 	      /* set cut rank */
1943 	      SCIProwChgRank(aggrcut, aggrcutrank); /*lint !e644*/
1944 	
1945 	      /* cache the row extension and only flush them if the cut gets added */
1946 	      SCIP_CALL( SCIPcacheRowExtensions(scip, aggrcut) );
1947 	
1948 	      /* collect all non-zero coefficients */
1949 	      for( int i = 0; i < aggrcutnnz; i++ )
1950 	      {
1951 	         var = SCIPcolGetVar(cols[aggrcutinds[i]]);
1952 	         SCIP_CALL( SCIPaddVarToRow(scip, aggrcut, var, aggrcutcoefs[i]) );
1953 	      }
1954 	
1955 	      /* flush all changes before adding the cut */
1956 	      SCIP_CALL( SCIPflushRowExtensions(scip, aggrcut) );
1957 	
1958 	      if( SCIProwGetNNonz(aggrcut) == 0 )
1959 	      {
1960 	         assert( SCIPisFeasNegative(scip, aggrcutrhs) );
1961 	         SCIPdebugMsg(scip, " -> Lagromory cut detected node infeasibility with cut 0 <= %g.\n", aggrcutrhs);
1962 	         *cutoff = TRUE;
1963 	         goto TERMINATE;
1964 	      }
1965 	      else
1966 	      {
1967 	         /* gathering lhs and rhs again in case the separator is extended later to add other cuts/constraints that may
1968 	          * have non-inf lhs or inf rhs */
1969 	         cutlhs = SCIProwGetLhs(aggrcut);
1970 	         cutrhs = SCIProwGetRhs(aggrcut);
1971 	         assert(SCIPisInfinity(scip, -cutlhs));
1972 	         assert(!SCIPisInfinity(scip, cutrhs));
1973 	
1974 	         SCIPdebugMsg(scip, " -> %s cut <%s>: rhs=%f\n", "lagromory", aggrcutname, aggrcutrhs);
1975 	
1976 	         /* check if the cut leads to infeasibility (i.e., *cutoff = TRUE) */
1977 	         {
1978 	            /* modifiable cuts cannot be declared infeasible, since we don't know all coefficients */
1979 	            if( SCIProwIsModifiable(aggrcut) )
1980 	               *cutoff = FALSE;
1981 	
1982 	            /* check for activity infeasibility */
1983 	            minactivity = SCIPgetRowMinActivity(scip, aggrcut);
1984 	            maxactivity = SCIPgetRowMaxActivity(scip, aggrcut);
1985 	
1986 	            if( (!SCIPisInfinity(scip, cutrhs) && SCIPisFeasPositive(scip, minactivity - cutrhs)) ||
1987 	                  (!SCIPisInfinity(scip, -cutlhs) && SCIPisFeasNegative(scip, maxactivity - cutlhs)) )
1988 	            {
1989 	               SCIPdebugMsg(scip, "cut <%s> is infeasible (sides=[%g,%g], act=[%g,%g])\n",
1990 	                     SCIProwGetName(aggrcut), cutlhs, cutrhs, minactivity, maxactivity);
1991 	               *cutoff = TRUE;
1992 	            }
1993 	         }
1994 	
1995 	         /* add the aggregated cut to a separate data structure */
1996 	         aggrcuts[*naggrcuts] = aggrcut;
1997 	         (*naggrcuts)++;
1998 	      }
1999 	
2000 	      QUAD_ASSIGN(quadterm, 0.0);
2001 	      for( int i = 0; i < ncols; i++ )
2002 	      {
2003 	         aggrcutvals[i] = 0.0;
2004 	         QUAD_ARRAY_STORE(tmpcutvals, i, quadterm);
2005 	      }
2006 	   }
2007 	
2008 	TERMINATE:
2009 	   /* free memory */
2010 	   SCIPfreeBufferArray(scip, &aggrcutinds);
2011 	   SCIPfreeBufferArray(scip, &aggrcutcoefs);
2012 	   SCIPfreeCleanBufferArray(scip, &aggrcutvals);
2013 	   SCIPfreeCleanBufferArray(scip, &tmpcutvals);
2014 	
2015 	   return SCIP_OKAY;
2016 	}
2017 	
2018 	/** main method: LP solution separation method of separator */
2019 	static
2020 	SCIP_RETCODE generateGMICuts(
2021 	   SCIP*                 scip,                     /**< SCIP data structure */
2022 	   SCIP_SEPA*            sepa,                     /**< pointer to the separator */
2023 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
2024 	   int                   mainiternum,              /**< iteration number of the outer loop of the relax-and-cut algorithm */
2025 	   int                   subgradientiternum,       /**< iteration number of the subgradient algorithm */
2026 	   SCIP_SOL*             sol,                      /**< LP solution to be used for cut generation */
2027 	   SCIP_Real*            solvals,                  /**< values of the LP solution to be used for cut generation */
2028 	   int                   nmaxgeneratedperroundcuts,/**< maximal number of cuts allowed to generate per separation round */
2029 	   SCIP_Bool             allowlocal,               /**< should locally valid cuts be generated? */
2030 	   SCIP_ROW**            generatedcurrroundcuts,   /**< cuts generated in the current separation round */
2031 	   SCIP_Real*            generatedcutefficacies,   /**< array of generated cut efficacies w.r.t. the respective LP bases used for cut
2032 	                                                     generations */
2033 	   int                   ngeneratedcurrroundcuts,  /**< number of cuts generated until the previous basis in the current separation round */
2034 	   int*                  ngeneratednewcuts,        /**< number of new cuts generated using the current basis */
2035 	   int                   depth,                    /**< depth of the current node in the tree */
2036 	   SCIP_Bool*            cutoff                    /**< should the current node be cutoff? */
2037 	   )
2038 	{
2039 	   SCIP_Real minfrac;
2040 	   SCIP_Real maxfrac;
2041 	   SCIP_Real frac;
2042 	   SCIP_Real* basisfrac;
2043 	   SCIP_Real* binvrow;
2044 	   SCIP_AGGRROW* aggrrow;
2045 	   SCIP_Bool success;
2046 	   SCIP_Real* cutcoefs;
2047 	   SCIP_Real cutrhs;
2048 	   SCIP_Real cutefficacy;
2049 	   SCIP_Bool cutislocal;
2050 	   SCIP_ROW** rows;
2051 	   SCIP_COL** cols;
2052 	   SCIP_VAR* var;
2053 	   SCIP_ROW* row;
2054 	   int cutrank;
2055 	   int cutnnz;
2056 	   int* cutinds;
2057 	   int* basisind;
2058 	   int* inds;
2059 	   int nrows;
2060 	   int ncols;
2061 	   int* basisperm;
2062 	   int k;
2063 	   int c;
2064 	   int ninds;
2065 	   int nmaxcutsperlp;
2066 	
2067 	   assert(ngeneratednewcuts != NULL);
2068 	
2069 	   minfrac = sepadata->away;
2070 	   maxfrac = 1.0 - sepadata->away;
2071 	   SCIP_CALL( SCIPgetLPRowsData(scip, &rows, &nrows) );
2072 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
2073 	   *ngeneratednewcuts = 0;
2074 	   nmaxcutsperlp = ((depth == 0) ? sepadata->nmaxcutsperlproot : sepadata->nmaxcutsperlp);
2075 	
2076 	   /* allocate memory */
2077 	   SCIP_CALL( SCIPallocBufferArray(scip, &basisperm, nrows) );
2078 	   SCIP_CALL( SCIPallocBufferArray(scip, &basisfrac, nrows) );
2079 	   SCIP_CALL( SCIPallocBufferArray(scip, &binvrow, nrows) );
2080 	   SCIP_CALL( SCIPallocBufferArray(scip, &inds, nrows) );
2081 	   SCIP_CALL( SCIPallocBufferArray(scip, &basisind, nrows) );
2082 	   SCIP_CALL( SCIPallocBufferArray(scip, &cutcoefs, ncols) );
2083 	   SCIP_CALL( SCIPallocBufferArray(scip, &cutinds, ncols) );
2084 	   SCIP_CALL( SCIPaggrRowCreate(scip, &aggrrow) );
2085 	
2086 	   SCIP_CALL( SCIPgetLPBasisInd(scip, basisind) );
2087 	
2088 	   /* check if the rows of the simplex tableau are suitable for cut generation and build an array of fractions */
2089 	   for( int i = 0; i < nrows; ++i )
2090 	   {
2091 	      frac = 0.0;
2092 	
2093 	      c = basisind[i];
2094 	
2095 	      basisperm[i] = i;
2096 	
2097 	      /* if the simplex tableau row corresponds to an LP column */
2098 	      if( c >= 0 )
2099 	      {
2100 	         assert(c < ncols);
2101 	
2102 	         var = SCIPcolGetVar(cols[c]);
2103 	         /* if the column is non-continuous one */
2104 	         if( SCIPvarGetType(var) != SCIP_VARTYPE_CONTINUOUS )
2105 	         {
2106 	            frac = SCIPfeasFrac(scip, solvals[c]);
2107 	            frac = MIN(frac, 1.0 - frac);
2108 	         }
2109 	      }
2110 	      /* if the simplex tableau row corresponds to an LP row and separation on rows is allowed */
2111 	      else if( sepadata->separaterows )
2112 	      {
2113 	         assert(0 <= -c-1);
2114 	         assert(-c-1 < nrows);
2115 	
2116 	         row = rows[-c-1];
2117 	         /* if the row is suitable for cut generation */
2118 	         if( SCIProwIsIntegral(row) && !SCIProwIsModifiable(row) )
2119 	         {
2120 	            frac = SCIPfeasFrac(scip, SCIPgetRowActivity(scip, row));
2121 	            frac = MIN(frac, 1.0 - frac);
2122 	         }
2123 	      }
2124 	
2125 	      if( frac >= minfrac )
2126 	      {
2127 	         /* slightly change fractionality to have random order for equal fractions */
2128 	         basisfrac[i] = frac + SCIPrandomGetReal(sepadata->randnumgen, -1e-6, 1e-6);
2129 	      }
2130 	      else
2131 	      {
2132 	         basisfrac[i] = 0.0;
2133 	      }
2134 	   }
2135 	
2136 	   /* if there is a need to sort the fractionalities */
2137 	   if( sepadata->sortcutoffsol )
2138 	   {
2139 	      /* sort basis indices by fractionality */
2140 	      SCIPsortDownRealInt(basisfrac, basisperm, nrows);
2141 	   }
2142 	
2143 	   /* for all basic columns belonging to integer variables, try to generate a GMI cut */
2144 	   for( int i = 0; i < nrows && !SCIPisStopped(scip) && !*cutoff; ++i )
2145 	   {
2146 	      if( (ngeneratedcurrroundcuts + *ngeneratednewcuts >= nmaxgeneratedperroundcuts) ||
2147 	            (sepadata->ntotalcuts + *ngeneratednewcuts >= sepadata->nmaxtotalcuts) ||
2148 	            (*ngeneratednewcuts >= nmaxcutsperlp) )
2149 	         break;
2150 	
2151 	      ninds = -1;
2152 	      cutefficacy = 0.0;
2153 	
2154 	      /* either break the loop or proceed to the next iteration if the fractionality is zero */
2155 	      if( basisfrac[i] == 0.0 )
2156 	      {
2157 	         if( sepadata->sortcutoffsol )
2158 	            break;
2159 	         else
2160 	            continue;
2161 	      }
2162 	
2163 	      k = basisperm[i];
2164 	
2165 	      /* get the row of B^-1 for this basic integer variable with fractional solution value and call aggregate function */
2166 	      SCIP_CALL( SCIPgetLPBInvRow(scip, k, binvrow, inds, &ninds) );
2167 	
2168 	      SCIP_CALL( SCIPaggrRowSumRows(scip, aggrrow, binvrow, inds, ninds, sepadata->sidetypebasis, allowlocal, 2,
2169 	               (int) MAXAGGRLEN(ncols), &success) );
2170 	
2171 	      if( !success )
2172 	         continue;
2173 	
2174 	      /* @todo Currently we are using the SCIPcalcMIR() function to compute the coefficients of the Gomory
2175 	       *       cut. Alternatively, we could use the direct version (see thesis of Achterberg formula (8.4)) which
2176 	       *       leads to cut a of the form \sum a_i x_i \geq 1. Rumor has it that these cuts are better.
2177 	       */
2178 	
2179 	      /* try to create GMI cut out of the aggregation row */
2180 	      SCIP_CALL( SCIPcalcMIR(scip, sol, POSTPROCESS, BOUNDSWITCH, USEVBDS, allowlocal, FIXINTEGRALRHS, NULL,
2181 	               NULL, minfrac, maxfrac, 1.0, aggrrow, cutcoefs, &cutrhs, cutinds, &cutnnz, &cutefficacy, &cutrank,
2182 	               &cutislocal, &success) );
2183 	
2184 	      if( success )
2185 	      {
2186 	         assert(allowlocal || !cutislocal); /*lint !e644*/
2187 	
2188 	         SCIP_CALL( constructCutRow(scip, sepa, sepadata, mainiternum, subgradientiternum, cutnnz, cutinds, cutcoefs,
2189 	                  cutefficacy, cutrhs, cutislocal, cutrank, generatedcurrroundcuts, generatedcutefficacies,
2190 	                  ngeneratedcurrroundcuts, ngeneratednewcuts, cutoff));
2191 	      }
2192 	   }
2193 	
2194 	   /* free temporary memory */
2195 	   SCIPfreeBufferArray(scip, &cutinds);
2196 	   SCIPfreeBufferArray(scip, &cutcoefs);
2197 	   SCIPfreeBufferArray(scip, &basisind);
2198 	   SCIPfreeBufferArray(scip, &inds);
2199 	   SCIPfreeBufferArray(scip, &binvrow);
2200 	   SCIPfreeBufferArray(scip, &basisfrac);
2201 	   SCIPfreeBufferArray(scip, &basisperm);
2202 	   SCIPaggrRowFree(scip, &aggrrow);
2203 	
2204 	   return SCIP_OKAY;
2205 	}
2206 	
2207 	/** update objective vector w.r.t. the fixed Lagrangian multipliers */
2208 	static
2209 	SCIP_RETCODE updateObjectiveVector(
2210 	   SCIP*                 scip,               /**< SCIP data structure */
2211 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers vector */
2212 	   SCIP_ROW**            cuts,               /**< cuts generated so far in the current separation round */
2213 	   int                   ncuts,              /**< number of cuts generated so far in the current separation round */
2214 	   SCIP_Real*            origobjcoefs,       /**< original objective function coefficients of the node linear relaxation */
2215 	   SCIP_Bool*            objvecsdiffer       /**< whether the updated objective function coefficients differ from the old ones */
2216 	   )
2217 	{
2218 	   SCIP_Real* objvals;
2219 	   SCIP_Real* prod;
2220 	   SCIP_Real* oldobjvals;
2221 	   SCIP_Real* cutvals;
2222 	   SCIP_COL** cutcols;
2223 	   SCIP_COL** cols;
2224 	   SCIP_VAR* var;
2225 	   int cutnnz;
2226 	   int collppos;
2227 	   int ncols;
2228 	
2229 	   assert(objvecsdiffer != NULL);
2230 	   assert(ncuts > 0);
2231 	
2232 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
2233 	
2234 	   /* allocate memory */
2235 	   SCIP_CALL( SCIPallocBufferArray(scip, &objvals, ncols) );
2236 	   SCIP_CALL( SCIPallocBufferArray(scip, &oldobjvals, ncols) );
2237 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &prod, ncols) );
2238 	   *objvecsdiffer = FALSE;
2239 	
2240 	   /* find the product of Lagrangian multipliers and cut coefficients */
2241 	   for( int i = 0; i < ncuts; i++ )
2242 	   {
2243 	      if( !SCIPisZero(scip, dualvector[i]) )
2244 	      {
2245 	         assert(!(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])) && SCIPisInfinity(scip, SCIProwGetRhs(cuts[i]))));
2246 	         assert(SCIPisInfinity(scip, -SCIProwGetLhs(cuts[i])));
2247 	         assert(!SCIPisInfinity(scip, SCIProwGetRhs(cuts[i])));
2248 	
2249 	         cutnnz = SCIProwGetNNonz(cuts[i]);
2250 	         assert(cutnnz <= ncols);
2251 	         cutvals = SCIProwGetVals(cuts[i]);
2252 	         cutcols = SCIProwGetCols(cuts[i]);
2253 	
2254 	         for( int j = 0; j < cutnnz; ++j )
2255 	         {
2256 	            collppos = SCIPcolGetLPPos(cutcols[j]);
2257 	            assert(collppos >= 0);
2258 	            assert(collppos <= ncols);
2259 	
2260 	            prod[collppos] += dualvector[i] * cutvals[j];
2261 	         }
2262 	      }
2263 	   }
2264 	
2265 	   /* change objective coefficients */
2266 	   for( int i = 0; i < ncols; i++ )
2267 	   {
2268 	      var = SCIPcolGetVar(cols[i]);
2269 	      oldobjvals[i] = SCIPgetVarObjDive(scip, var);
2270 	      objvals[i] = origobjcoefs[i] + prod[i];
2271 	
2272 	      SCIP_CALL( SCIPchgVarObjDive(scip, var, objvals[i]) );
2273 	
2274 	      /* identify if the updated objective vector is indeed different from the previous one */
2275 	      if( !(*objvecsdiffer) && !SCIPisEQ(scip, oldobjvals[i], objvals[i]) )
2276 	         *objvecsdiffer = TRUE;
2277 	   }
2278 	
2279 	   for( int i = 0; i < ncols; i++)
2280 	   {
2281 	      prod[i] = 0.0;
2282 	   }
2283 	
2284 	   /* free memory  */
2285 	   SCIPfreeCleanBufferArray(scip, &prod);
2286 	   SCIPfreeBufferArray(scip, &oldobjvals);
2287 	   SCIPfreeBufferArray(scip, &objvals);
2288 	
2289 	   return SCIP_OKAY;
2290 	}
2291 	
2292 	/** add GMI cuts to the objective function of the Lagrangian dual problem by introducing new Lagrangian multipliers */
2293 	static
2294 	SCIP_RETCODE addGMICutsAsSoftConss(
2295 	   SCIP_Real*            dualvector,         /**< Lagrangian multipliers vector */
2296 	   int                   ngeneratedcuts,     /**< number of cuts generated so far in the current separation round */
2297 	   int*                  naddedcuts,         /**< number of cuts added so far in the current separation round to the Lagrangian dual problem
2298 	                                               upon penalization */
2299 	   int*                  nnewaddedsoftcuts   /**< number of cuts added newly to the Lagrangian dual problem upon penalization */
2300 	   )
2301 	{
2302 	   assert(*naddedcuts <= ngeneratedcuts);
2303 	
2304 	   /* set the initial penalty of the newly penalized cuts as zero */
2305 	   for( int i = *naddedcuts; i < ngeneratedcuts; i++ )
2306 	      dualvector[i] = 0.0;
2307 	
2308 	   *nnewaddedsoftcuts = ngeneratedcuts - *naddedcuts;
2309 	   *naddedcuts = ngeneratedcuts;
2310 	
2311 	   return SCIP_OKAY;
2312 	}
2313 	
2314 	/** solve the Lagrangian dual problem */
2315 	static
2316 	SCIP_RETCODE solveLagrangianDual(
2317 	   SCIP*                 scip,                     /**< SCIP data structure */
2318 	   SCIP_SEPA*            sepa,                     /**< pointer to the separator */
2319 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
2320 	   SCIP_SOL*             sol,                      /**< data structure to store an LP solution upon solving a Lagrangian dual problem with
2321 	                                                     fixed Lagrangian multipliers */
2322 	   SCIP_Real*            solvals,                  /**< values of the LP solution obtained upon solving a Lagrangian dual problem with fixed
2323 	                                                     Lagrangian multipliers */
2324 	   int                   mainiternum,              /**< iteration number of the outer loop of the relax-and-cut algorithm */
2325 	   SCIP_Real             ubparam,                  /**< estimate of the optimal Lagrangian dual value */
2326 	   int                   depth,                    /**< depth of the current node in the tree */
2327 	   SCIP_Bool             allowlocal,               /**< should locally valid cuts be generated? */
2328 	   int                   nmaxgeneratedperroundcuts,/**< maximal number of cuts allowed to generate per separation round */
2329 	   SCIP_Real*            origobjcoefs,             /**< original objective function coefficients of the node linear relaxation */
2330 	   SCIP_Real             origobjoffset,            /**< original objective function offset of the node linear relaxation */
2331 	   SCIP_Real*            dualvector,               /**< Lagrangian multipliers vector */
2332 	   int*                  nsoftcuts,                /**< number of generated cuts that were penalized and added to the Lagrangian dual problem */
2333 	   SCIP_ROW**            generatedcurrroundcuts,   /**< cuts generated in the current separation round */
2334 	   SCIP_Real*            generatedcutefficacies,   /**< array of generated cut efficacies w.r.t. the respective LP bases used for cut
2335 	                                                     generations */
2336 	   int*                  ngeneratedcutsperiter,    /**< number of cuts generated per subgradient iteration in the current separation round */
2337 	   int*                  ngeneratedcurrroundcuts,  /**< number of cuts generated so far in the current separation round */
2338 	   int*                  ncurrroundlpiters,        /**< number of LP iterations taken for solving Lagrangian dual problems with fixed
2339 	                                                     multipliers in the current separator round */
2340 	   SCIP_Bool*            cutoff,                   /**< should the current node be cutoff? */
2341 	   SCIP_Real*            bestlagrangianval,        /**< best Lagrangian value found so far */
2342 	   SCIP_Real*            bestdualvector,           /**< Lagrangian multipliers corresponding to the best Lagrangian value found so far */
2343 	   int*                  bestdualvectorlen,        /**< length of the Lagrangian multipliers corresponding to the best Lagrangian value
2344 	                                                     found so far */
2345 	   int*                  nbestdualupdates,         /**< number of best Lagrangian values found so far */
2346 	   int*                  totaliternum              /**< total number of iterations of the relax-and-cut algorithm performed so far */
2347 	   )
2348 	{
2349 	   SCIP_Real* subgradient;
2350 	   SCIP_Real muparam;
2351 	   SCIP_Real steplength;
2352 	   SCIP_Real objval;
2353 	   SCIP_Real lagrangianval;
2354 	   SCIP_Real* lagrangianvals;
2355 	   SCIP_Real avglagrangianval;
2356 	   SCIP_Real maxsoftcutviol;
2357 	   SCIP_Real maxnzsubgradientdualprod;
2358 	   SCIP_Real maxviolscore;
2359 	   SCIP_Real maxviolscoreold;
2360 	   SCIP_Real nviolscore;
2361 	   SCIP_Real nviolscoreold;
2362 	   SCIP_Real scoreweight;
2363 	   SCIP_Real ballradius;
2364 	   SCIP_Bool solfound;
2365 	   SCIP_Bool backtrack;
2366 	   SCIP_Bool terminate;
2367 	   SCIP_Bool subgradientzero;
2368 	   SCIP_Bool objvecsdiffer;
2369 	   SCIP_Bool dualvecsdiffer;
2370 	   SCIP_Bool solvelp;
2371 	   int ncurrroundlpiterslast;
2372 	   int nlpiters;
2373 	   int ngeneratednewcuts;
2374 	   int nnewaddedsoftcuts;
2375 	   int nsoftcutviols;
2376 	   int nnzsubgradientdualprod;
2377 	
2378 	   SCIP_CALL( SCIPallocBufferArray(scip, &lagrangianvals, sepadata->nmaxsubgradientiters) );
2379 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &subgradient, nmaxgeneratedperroundcuts) );
2380 	
2381 	   muparam = sepadata->muparaminit;
2382 	   steplength = 0.0;
2383 	   objval = 0.0;
2384 	   avglagrangianval = 0.0;
2385 	   maxsoftcutviol = 0.0;
2386 	   maxnzsubgradientdualprod = 0.0;
2387 	   maxviolscore = 0.0;
2388 	   nviolscore = 0.0;
2389 	   scoreweight = 1.0;
2390 	   ballradius = sepadata->radiusinit;
2391 	   ngeneratednewcuts = 0;
2392 	   nsoftcutviols = 0;
2393 	   nnzsubgradientdualprod = 0;
2394 	   terminate = FALSE;
2395 	   subgradientzero = FALSE;
2396 	   objvecsdiffer = FALSE;
2397 	   dualvecsdiffer = FALSE;
2398 	   solvelp = TRUE;
2399 	
2400 	   /* update objective vector based on input Lagrangian multipliers */
2401 	   if( *nsoftcuts > 0 )
2402 	   {
2403 	      SCIP_CALL( updateObjectiveVector(scip, dualvector, generatedcurrroundcuts, *nsoftcuts, origobjcoefs, &objvecsdiffer) );
2404 	   }
2405 	
2406 	   /* termination check */
2407 	   SCIP_CALL( checkLagrangianDualTermination(sepadata, -1, -1, FALSE, *ngeneratedcurrroundcuts,
2408 	            nmaxgeneratedperroundcuts, *ncurrroundlpiters, depth, &terminate) );
2409 	
2410 	   /* the subgradient algorithm loop */
2411 	   for( int i = 0; i < sepadata->nmaxsubgradientiters && !SCIPisStopped(scip) && !terminate; i++ )
2412 	   {
2413 	      solfound = FALSE;
2414 	      subgradientzero = FALSE;
2415 	      objvecsdiffer = FALSE;
2416 	      dualvecsdiffer = FALSE;
2417 	      nnewaddedsoftcuts = 0;
2418 	      scoreweight *= sepadata->radiusupdateweight;
2419 	
2420 	      ncurrroundlpiterslast = *ncurrroundlpiters;
2421 	      if( solvelp )
2422 	      {
2423 	         /* solve Lagrangian dual for fixed Lagrangian multipliers */
2424 	         SCIP_CALL( solveLagromoryLP(scip, sepadata, depth, origobjoffset, &solfound, sol, solvals, &objval,
2425 	                  ncurrroundlpiters) );
2426 	      }
2427 	      nlpiters = *ncurrroundlpiters - ncurrroundlpiterslast;
2428 	
2429 	      /* if an optimal solution has been found, then generate cuts and do other operations */
2430 	      if( solfound )
2431 	      {
2432 	         /* generate GMI cuts if a new basis solution is found */
2433 	         if( (nlpiters >= 1) && (i % sepadata->cutgenfreq == 0) )
2434 	         {
2435 	            ngeneratednewcuts = 0;
2436 	            SCIP_CALL( generateGMICuts(scip, sepa, sepadata, mainiternum, i, sol, solvals,
2437 	                     nmaxgeneratedperroundcuts, allowlocal, generatedcurrroundcuts, generatedcutefficacies,
2438 	                     *ngeneratedcurrroundcuts, &ngeneratednewcuts, depth, cutoff));
2439 	            sepadata->ntotalcuts += ngeneratednewcuts;
2440 	            *ngeneratedcurrroundcuts += ngeneratednewcuts;
2441 	            ngeneratedcutsperiter[mainiternum * sepadata->nmaxsubgradientiters + i + 1] = ngeneratednewcuts;
2442 	         }
2443 	
2444 	         /* update subgradient, i.e., find the residuals of the penalized cuts, and determine various violations */
2445 	         updateSubgradient(scip, sol, generatedcurrroundcuts, *nsoftcuts, subgradient, dualvector, &subgradientzero,
2446 	               &nsoftcutviols, &maxsoftcutviol, &nnzsubgradientdualprod, &maxnzsubgradientdualprod);
2447 	
2448 	         /* calculate Lagrangian value for the fixed Lagrangian multipliers, and update best and avg values */
2449 	         updateLagrangianValue(scip, objval, dualvector, generatedcurrroundcuts, *nsoftcuts, &lagrangianval);
2450 	         if( SCIPisPositive(scip, lagrangianval - *bestlagrangianval) )
2451 	         {
2452 	            *bestlagrangianval = lagrangianval;
2453 	            for( int j = 0; j < *nsoftcuts; j++ )
2454 	            {
2455 	               bestdualvector[j] = dualvector[j];
2456 	            }
2457 	            *bestdualvectorlen = *nsoftcuts;
2458 	            (*nbestdualupdates)++;
2459 	         }
2460 	         lagrangianvals[i] = lagrangianval;
2461 	         if( i < sepadata->nmaxlagrangianvalsforavg )
2462 	         {
2463 	            avglagrangianval = (avglagrangianval * i + lagrangianval)/(i+1);
2464 	         }
2465 	         else
2466 	         {
2467 	            avglagrangianval = (avglagrangianval * sepadata->nmaxlagrangianvalsforavg -
2468 	                  lagrangianvals[i - sepadata->nmaxlagrangianvalsforavg] +
2469 	                  lagrangianval)/(sepadata->nmaxlagrangianvalsforavg);
2470 	         }
2471 	
2472 	         /* if the subgradient vector is non-zero, then update the mu parameter and the Lagrangian multipliers */
2473 	         if( !subgradientzero )
2474 	         {
2475 	            /* update mu param */
2476 	            SCIP_CALL( updateMuSteplengthParam(scip, sepadata, i, ubparam, lagrangianvals, *bestlagrangianval, avglagrangianval,
2477 	                     &muparam, &backtrack) );
2478 	
2479 	            /* update step length */
2480 	            updateStepLength(scip, muparam, ubparam, lagrangianval, subgradient, *nsoftcuts, &steplength);
2481 	
2482 	            /* update scores to determine how to update the stabilization ball radius */
2483 	            maxviolscoreold = maxviolscore;
2484 	            nviolscoreold = nviolscore;
2485 	            maxviolscore = (1.0 - scoreweight) * maxsoftcutviol + scoreweight * maxnzsubgradientdualprod;
2486 	            nviolscore = (1.0 - scoreweight) * nsoftcutviols + scoreweight * nnzsubgradientdualprod;
2487 	
2488 	            /* update Lagrangian multipliers */
2489 	            SCIP_CALL( updateDualVector(scip, sepadata, dualvector, bestdualvector, *bestdualvectorlen,
2490 	                     *nbestdualupdates, i, *totaliternum, steplength, subgradient, *nsoftcuts, backtrack, maxviolscore,
2491 	                     maxviolscoreold, nviolscore, nviolscoreold, nlpiters, &dualvecsdiffer, &ballradius) );
2492 	
2493 	            /* update objective vector based on updated Lagrangian multipliers */
2494 	            if( dualvecsdiffer )
2495 	            {
2496 	               SCIP_CALL( updateObjectiveVector(scip, dualvector, generatedcurrroundcuts, *nsoftcuts, origobjcoefs, &objvecsdiffer) );
2497 	            }
2498 	         }
2499 	         /* if the subgradient vector if zero, then simply mark that the Lagrangian multipliers and the objective
2500 	          * function of the Lagrangian dual did not change */
2501 	         else
2502 	         {
2503 	            dualvecsdiffer = FALSE;
2504 	            objvecsdiffer = FALSE;
2505 	         }
2506 	
2507 	
2508 	         /* add generated GMI cuts to the objective function of the Lagrangian dual problem by introducing new
2509 	          * Lagrangian multipliers */
2510 	         if( (i % sepadata->cutaddfreq == 0) || (!dualvecsdiffer && !objvecsdiffer &&
2511 	                  (*ngeneratedcurrroundcuts - *nsoftcuts > 0)) )
2512 	         {
2513 	            SCIP_CALL( addGMICutsAsSoftConss(dualvector, *ngeneratedcurrroundcuts, nsoftcuts, &nnewaddedsoftcuts) );
2514 	         }
2515 	      }
2516 	      else
2517 	      {
2518 	         /* add any remaining generated GMI cuts to the objective function of the Lagrangian dual problem by introducing
2519 	          * new Lagrangian multipliers */
2520 	         if( (*ngeneratedcurrroundcuts - *nsoftcuts) > 0 )
2521 	         {
2522 	            SCIP_CALL( addGMICutsAsSoftConss(dualvector, *ngeneratedcurrroundcuts, nsoftcuts, &nnewaddedsoftcuts) );
2523 	         }
2524 	
2525 	         solvelp = FALSE;
2526 	      }
2527 	
2528 	      /* termination check */
2529 	      SCIP_CALL( checkLagrangianDualTermination(sepadata, nnewaddedsoftcuts, *ngeneratedcurrroundcuts - *nsoftcuts,
2530 	               objvecsdiffer, *ngeneratedcurrroundcuts, nmaxgeneratedperroundcuts, *ncurrroundlpiters, depth,
2531 	               &terminate) );
2532 	
2533 	      (*totaliternum)++;
2534 	   }
2535 	
2536 	   /* add any remaining generated GMI cuts to the objective function of the Lagrangian dual problem by introducing new
2537 	    * Lagrangian multipliers */
2538 	   if( (*ngeneratedcurrroundcuts - *nsoftcuts) > 0 )
2539 	   {
2540 	      SCIP_CALL( addGMICutsAsSoftConss(dualvector, *ngeneratedcurrroundcuts, nsoftcuts, &nnewaddedsoftcuts) );
2541 	   }
2542 	
2543 	   /* free memory */
2544 	   for( int i = 0; i < nmaxgeneratedperroundcuts; i++)
2545 	   {
2546 	      subgradient[i] = 0.0;
2547 	   }
2548 	   SCIPfreeCleanBufferArray(scip, &subgradient);
2549 	   SCIPfreeBufferArray(scip, &lagrangianvals);
2550 	
2551 	   return SCIP_OKAY;
2552 	}
2553 	
2554 	/** generates initial cut pool before solving the Lagrangian dual */
2555 	static
2556 	SCIP_RETCODE generateInitCutPool(
2557 	   SCIP*                 scip,                     /**< SCIP data structure */
2558 	   SCIP_SEPA*            sepa,                     /**< separator */
2559 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
2560 	   int                   mainiternum,              /**< iteration number of the outer loop of the relax-and-cut algorithm */
2561 	   SCIP_SOL*             sol,                      /**< LP solution to be used for cut generation */
2562 	   SCIP_Real*            solvals,                  /**< values of the LP solution to be used for cut generation */
2563 	   int                   nmaxgeneratedperroundcuts,/**< maximal number of cuts allowed to generate per separation round */
2564 	   SCIP_Bool             allowlocal,               /**< should locally valid cuts be generated? */
2565 	   SCIP_ROW**            generatedcurrroundcuts,   /**< cuts generated in the current separation round */
2566 	   SCIP_Real*            generatedcutefficacies,   /**< array of generated cut efficacies w.r.t. the respective LP bases used for cut
2567 	                                                     generations */
2568 	   int*                  ngeneratedcutsperiter,    /**< number of cuts generated per subgradient iteration in the current separation round */
2569 	   int*                  ngeneratedcurrroundcuts,  /**< number of cuts generated so far in the current separation round */
2570 	   int                   depth,                    /**< depth of the current node in the tree */
2571 	   SCIP_Bool*            cutoff                    /**< should the current node be cutoff? */
2572 	   )
2573 	{
2574 	   int ngeneratednewcuts;
2575 	
2576 	   ngeneratednewcuts = 0;
2577 	
2578 	   /* generate initial set of cuts */
2579 	   SCIP_CALL( generateGMICuts(scip, sepa, sepadata, mainiternum, -1, sol, solvals, nmaxgeneratedperroundcuts,
2580 	            allowlocal, generatedcurrroundcuts, generatedcutefficacies, *ngeneratedcurrroundcuts, &ngeneratednewcuts,
2581 	            depth, cutoff) );
2582 	
2583 	   /* update certain statistics */
2584 	   sepadata->ntotalcuts += ngeneratednewcuts;
2585 	   *ngeneratedcurrroundcuts += ngeneratednewcuts;
2586 	   ngeneratedcutsperiter[sepadata->nmaxsubgradientiters * mainiternum] = ngeneratednewcuts;
2587 	
2588 	   return SCIP_OKAY;
2589 	}
2590 	
2591 	/** add cuts to SCIP */
2592 	static
2593 	SCIP_RETCODE addCuts(
2594 	   SCIP*                 scip,               /**< SCIP data structure */
2595 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
2596 	   SCIP_ROW**            cuts,               /**< cuts generated so far in the current separation round */
2597 	   int                   ncuts,              /**< number of cuts generated so far in the current separation round */
2598 	   SCIP_Longint          maxdnom,            /**< maximum denominator in the rational representation of cuts */
2599 	   SCIP_Real             maxscale,           /**< maximal scale factor to scale the cuts to integral values */
2600 	   int*                  naddedcuts,         /**< number of cuts added to either global cutpool or sepastore */
2601 	   SCIP_Bool*            cutoff              /**< should the current node be cutoff? */
2602 	   )
2603 	{
2604 	   SCIP_ROW* cut;
2605 	   SCIP_Bool madeintegral;
2606 	
2607 	   assert(scip != NULL);
2608 	   assert(sepadata != NULL);
2609 	   assert(naddedcuts != NULL);
2610 	   assert(cutoff != NULL);
2611 	
2612 	   *cutoff = FALSE;
2613 	   madeintegral = FALSE;
2614 	
2615 	   for( int i = 0; i < ncuts && !*cutoff; i++ )
2616 	   {
2617 	      cut = cuts[i];
2618 	
2619 	      if( SCIProwGetNNonz(cut) == 1 )
2620 	      {
2621 	         /* Add the bound change as cut to avoid that the LP gets modified. This would mean that the LP is not flushed
2622 	          * and the method SCIPgetLPBInvRow() fails; SCIP internally will apply this bound change automatically. */
2623 	         SCIP_CALL( SCIPaddRow(scip, cut, TRUE, cutoff) );
2624 	         ++(*naddedcuts);
2625 	      }
2626 	      else
2627 	      {
2628 	         if( sepadata->makeintegral && SCIPgetRowNumIntCols(scip, cut) == SCIProwGetNNonz(cut) )
2629 	         {
2630 	            /* try to scale the cut to integral values */
2631 	            SCIP_CALL( SCIPmakeRowIntegral(scip, cut, -SCIPepsilon(scip), SCIPsumepsilon(scip),
2632 	                     maxdnom, maxscale, MAKECONTINTEGRAL, &madeintegral) );
2633 	
2634 	            /* if RHS = plus infinity (due to scaling), the cut is useless, so we are not adding it */
2635 	            if( madeintegral && SCIPisInfinity(scip, SCIProwGetRhs(cut)) )
2636 	               return SCIP_OKAY;
2637 	         }
2638 	
2639 	         if( SCIPisCutNew(scip, cut) )
2640 	         {
2641 	            /* add global cuts which are not implicit bound changes to the cut pool */
2642 	            if( !SCIProwIsLocal(cut) )
2643 	            {
2644 	               if( sepadata->delayedcuts )
2645 	               {
2646 	                  SCIP_CALL( SCIPaddDelayedPoolCut(scip, cut) );
2647 	               }
2648 	               else
2649 	               {
2650 	                  SCIP_CALL( SCIPaddPoolCut(scip, cut) );
2651 	               }
2652 	            }
2653 	            else
2654 	            {
2655 	               /* local cuts we add to the sepastore */
2656 	               SCIP_CALL( SCIPaddRow(scip, cut, sepadata->forcecuts, cutoff) );
2657 	            }
2658 	            ++(*naddedcuts);
2659 	         }
2660 	      }
2661 	   }
2662 	
2663 	   return SCIP_OKAY;
2664 	}
2665 	
2666 	/** check different termination criteria */
2667 	/* @todo nlpssolved criterion? */
2668 	static
2669 	SCIP_RETCODE checkMainLoopTermination(
2670 	   SCIP_SEPADATA*        sepadata,                 /**< separator data structure */
2671 	   SCIP_Bool             cutoff,                   /**< should the current node be cutoff? */
2672 	   SCIP_Bool             dualvecsdiffer,           /**< whether the updated Lagrangian multipliers differ from the old one */
2673 	   int                   ngeneratedcurrroundcuts,  /**< number of cuts generated in the current separation round */
2674 	   int                   nsoftcuts,                /**< number of generated cuts that were penalized and added to the Lagrangian dual problem */
2675 	   int                   nmaxgeneratedperroundcuts,/**< maximal number of cuts allowed to generate per separation round */
2676 	   int                   ncurrroundlpiters,        /**< number of LP iterations taken for solving Lagrangian dual problems with fixed
2677 	                                                     multipliers in the current separator round */
2678 	   int                   depth,                    /**< depth of the current node in the tree */
2679 	   SCIP_Bool*            terminate                 /**< whether to terminate the relax-and-cut algorithm */
2680 	   )
2681 	{
2682 	   *terminate = FALSE;
2683 	
2684 	   /* check if the node has been identified to be cutoff */
2685 	   if( cutoff )
2686 	      *terminate = TRUE;
2687 	
2688 	   /* check if the Lagrangian multipliers do not differ from the previous iteration and no new cuts exist for penalizing */
2689 	   if( !dualvecsdiffer && (ngeneratedcurrroundcuts == nsoftcuts) )
2690 	      *terminate = TRUE;
2691 	
2692 	   /* check if allowed number of cuts in this separation round have already been generated */
2693 	   if( ngeneratedcurrroundcuts >= nmaxgeneratedperroundcuts )
2694 	      *terminate = TRUE;
2695 	
2696 	   /* check if allowed number of cuts in the tree have already been generated */
2697 	   if( sepadata->ntotalcuts >= sepadata->nmaxtotalcuts )
2698 	      *terminate = TRUE;
2699 	
2700 	   /* check if allowed number of simplex iterations in this separation round have already been used up */
2701 	   if( (sepadata->nmaxperroundlpiters >= 0) && (ncurrroundlpiters >= sepadata->nmaxperroundlpiters) )
2702 	      *terminate = TRUE;
2703 	
2704 	   /* check if allowed number of simplex iterations in the root node have already been used up */
2705 	   if( (depth == 0) && (sepadata->nmaxrootlpiters >= 0) && (sepadata->nrootlpiters >= sepadata->nmaxrootlpiters) )
2706 	      *terminate = TRUE;
2707 	
2708 	   /* check if allowed number of simplex iterations in the tree have already been used up */
2709 	   if( (sepadata->nmaxtotallpiters >= 0) && (sepadata->ntotallpiters >= sepadata->nmaxtotallpiters) )
2710 	      *terminate = TRUE;
2711 	
2712 	   return SCIP_OKAY;
2713 	}
2714 	
2715 	/** Searches and tries to add Lagromory cuts */
2716 	static
2717 	SCIP_RETCODE separateCuts(
2718 	   SCIP*                 scip,               /**< SCIP data structure */
2719 	   SCIP_SEPA*            sepa,               /**< separator */
2720 	   SCIP_SEPADATA*        sepadata,           /**< separator data structure */
2721 	   SCIP_Real             ubparam,            /**< estimate of the optimal Lagrangian dual value */
2722 	   int                   depth,              /**< depth of the current node in the tree */
2723 	   SCIP_Bool             allowlocal,         /**< should locally valid cuts be generated? */
2724 	   SCIP_RESULT*          result              /**< final result of the separation round */
2725 	   )
2726 	{
2727 	   SCIP_ROW** generatedcurrroundcuts;
2728 	   SCIP_ROW** aggregatedcurrroundcuts;
2729 	   SCIP_Real* generatedcutefficacies;
2730 	   SCIP_Bool solfound;
2731 	   SCIP_SOL* softcutslpsol;
2732 	   SCIP_Real* softcutslpsolvals;
2733 	   SCIP_SOL* hardcutslpsol;
2734 	   SCIP_Real* hardcutslpsolvals;
2735 	   SCIP_Real* dualsol;
2736 	   SCIP_Real* dualvector;
2737 	   SCIP_Real* bestdualvector;
2738 	   SCIP_Real bestlagrangianval;
2739 	   SCIP_Real* origobjcoefs;
2740 	   SCIP_Real origobjoffset;
2741 	   SCIP_Real objval;
2742 	   SCIP_Real maxscale;
2743 	   SCIP_Longint maxdnom;
2744 	   SCIP_Bool cutoff;
2745 	   SCIP_Bool cutoff2;
2746 	   SCIP_Bool dualvecsdiffer;
2747 	   SCIP_Bool terminate;
2748 	   SCIP_COL** cols;
2749 	
2750 	   int* ngeneratedcutsperiter;
2751 	   int bestdualvectorlen;
2752 	   int nbestdualupdates;
2753 	   int totaliternum;
2754 	   int* cutindsperm;
2755 	   int nprocessedcuts;
2756 	   int ncols;
2757 	   int nrows;
2758 	   int ngeneratedcurrroundcuts;
2759 	   int nselectedcurrroundcuts;
2760 	   int nselectedcuts;
2761 	   int naddedcurrroundcuts;
2762 	   int naggregatedcurrroundcuts;
2763 	   int nmaxgeneratedperroundcuts;
2764 	   int ncurrroundlpiters;
2765 	   int nsoftcuts;
2766 	   int nsoftcutsold;
2767 	   int maxdepth;
2768 	
2769 	   assert(*result == SCIP_DIDNOTRUN);
2770 	   assert(sepadata != NULL);
2771 	   sepadata->ncalls = SCIPsepaGetNCalls(sepa);
2772 	   objval = SCIPinfinity(scip);
2773 	   bestlagrangianval = -SCIPinfinity(scip);
2774 	   bestdualvectorlen = 0;
2775 	   nbestdualupdates = 0;
2776 	   totaliternum = 0;
2777 	   ngeneratedcurrroundcuts = 0;
2778 	   naddedcurrroundcuts = 0;
2779 	   naggregatedcurrroundcuts = 0;
2780 	   ncurrroundlpiters = 0;
2781 	   nsoftcuts = 0;
2782 	   solfound = FALSE;
2783 	   cutoff = FALSE;
2784 	   cutoff2 = FALSE;
2785 	   dualvecsdiffer = FALSE;
2786 	   terminate = FALSE;
2787 	
2788 	   SCIPdebugMsg(scip, "Separating cuts...\n");
2789 	
2790 	   /* initialize the LP that will be used to solve the Lagrangian dual with fixed Lagrangian multipliers */
2791 	   SCIP_CALL( createLPWithSoftCuts(scip, sepadata) );
2792 	   /* create the LP that represents the node relaxation including all the generated cuts in this separator */
2793 	   SCIP_CALL( createLPWithHardCuts(scip, sepadata, NULL, 0) );
2794 	
2795 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
2796 	   nrows = SCIPgetNLPRows(scip);
2797 	
2798 	   /* get the maximal number of cuts allowed in a separation round */
2799 	   nmaxgeneratedperroundcuts = ((depth == 0) ? sepadata->nmaxperroundcutsroot : sepadata->nmaxperroundcuts);
2800 	
2801 	   /* allocate memory */
2802 	   SCIP_CALL( SCIPallocBufferArray(scip, &generatedcurrroundcuts, nmaxgeneratedperroundcuts) );
2803 	   SCIP_CALL( SCIPallocBufferArray(scip, &aggregatedcurrroundcuts, nmaxgeneratedperroundcuts) );
2804 	   SCIP_CALL( SCIPallocBufferArray(scip, &generatedcutefficacies, nmaxgeneratedperroundcuts) );
2805 	   SCIP_CALL( SCIPallocBufferArray(scip, &cutindsperm, nmaxgeneratedperroundcuts) );
2806 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &ngeneratedcutsperiter, sepadata->nmaxmainiters *
2807 	            (sepadata->nmaxsubgradientiters + 1)) );
2808 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &dualsol, nrows + nmaxgeneratedperroundcuts) );
2809 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &dualvector, nmaxgeneratedperroundcuts) );
2810 	   SCIP_CALL( SCIPallocCleanBufferArray(scip, &bestdualvector, nmaxgeneratedperroundcuts) );
2811 	   SCIP_CALL( SCIPallocBufferArray(scip, &softcutslpsolvals, ncols) );
2812 	   SCIP_CALL( SCIPcreateSol(scip, &softcutslpsol, NULL) );
2813 	   SCIP_CALL( SCIPallocBufferArray(scip, &hardcutslpsolvals, ncols) );
2814 	   SCIP_CALL( SCIPcreateSol(scip, &hardcutslpsol, NULL) );
2815 	   SCIP_CALL( SCIPallocBufferArray(scip, &origobjcoefs, ncols) );
2816 	
2817 	   /* store current objective function */
2818 	   for( int i = 0; i < ncols; i++ )
2819 	   {
2820 	      origobjcoefs[i] = SCIPcolGetObj(cols[i]);
2821 	   }
2822 	   origobjoffset = SCIPgetTransObjoffset(scip);
2823 	
2824 	   /* solve node LP relaxation to have an initial simplex tableau */
2825 	   SCIP_CALL( solveLagromoryLP(scip, sepadata, depth, origobjoffset, &solfound, softcutslpsol, softcutslpsolvals, &objval,
2826 	            &ncurrroundlpiters));
2827 	
2828 	   /* generate initial cut pool */
2829 	   SCIP_CALL( generateInitCutPool(scip, sepa, sepadata, 0, softcutslpsol, softcutslpsolvals, nmaxgeneratedperroundcuts, allowlocal,
2830 	            generatedcurrroundcuts, generatedcutefficacies, ngeneratedcutsperiter, &ngeneratedcurrroundcuts, depth,
2831 	            &cutoff) );
2832 	
2833 	   /* termination check */
2834 	   SCIP_CALL( checkMainLoopTermination(sepadata, cutoff, TRUE, ngeneratedcurrroundcuts, nsoftcuts,
2835 	            nmaxgeneratedperroundcuts, ncurrroundlpiters, depth, &terminate) );
2836 	
2837 	   /* compute cuts for each integer col with fractional val */
2838 	   for( int i = 0; i < sepadata->nmaxmainiters && !SCIPisStopped(scip) && !terminate; ++i )
2839 	   {
2840 	      nsoftcutsold = nsoftcuts;
2841 	      nsoftcuts = ngeneratedcurrroundcuts;
2842 	      dualvecsdiffer = FALSE;
2843 	
2844 	      /* solve Lagrangain dual */
2845 	      SCIP_CALL( solveLagrangianDual(scip, sepa, sepadata, softcutslpsol, softcutslpsolvals, i, ubparam, depth, allowlocal,
2846 	               nmaxgeneratedperroundcuts, origobjcoefs, origobjoffset, dualvector, &nsoftcuts, generatedcurrroundcuts,
2847 	               generatedcutefficacies, ngeneratedcutsperiter, &ngeneratedcurrroundcuts, &ncurrroundlpiters, &cutoff,
2848 	               &bestlagrangianval, bestdualvector, &bestdualvectorlen, &nbestdualupdates, &totaliternum) );
2849 	
2850 	      /* @todo filter cuts before adding them to the new LP that was created based on the node relaxation? */
2851 	
2852 	      /* update the LP representing the node relaxation by adding newly generated cuts */
2853 	      if( !cutoff && (ngeneratedcurrroundcuts - nsoftcutsold > 0) )
2854 	      {
2855 	         SCIP_CALL( createLPWithHardCuts(scip, sepadata, generatedcurrroundcuts, ngeneratedcurrroundcuts) );
2856 	
2857 	         /* solve the LP and get relevant information */
2858 	         SCIP_CALL( solveLPWithHardCuts(scip, sepadata, &solfound, hardcutslpsol, hardcutslpsolvals) );
2859 	
2860 	         /* if primal solution is found, then pass it to trysol heuristic */
2861 	         /* @note if trysol heuristic is was not present, then the solution found above, which can potentially be a MIP
2862 	          * primal feasible solution, will go to waste.
2863 	          */
2864 	         if( solfound && sepadata->heurtrysol != NULL )
2865 	         {
2866 	            SCIP_CALL( SCIPheurPassSolTrySol(scip, sepadata->heurtrysol, hardcutslpsol) );
2867 	         }
2868 	
2869 	         /* if dual feasible, then fetch dual solution and reset Lagrangian multipliers based on it. otherwise, retain the
2870 	          * Lagrangian multipliers and simply initialize the new multipliers to zeroes. */
2871 	         if( SCIPlpiIsDualFeasible(sepadata->lpiwithhardcuts) )
2872 	         {
2873 	            SCIP_CALL( SCIPlpiGetSol(sepadata->lpiwithhardcuts, NULL, NULL, dualsol, NULL, NULL) );
2874 	            SCIP_CALL( updateDualVector(scip, sepadata, dualvector, &(dualsol[nrows]),
2875 	                     ngeneratedcurrroundcuts, 0, -1, -1, 0.0, NULL, ngeneratedcurrroundcuts, TRUE, 0.0, 0.0, 0.0, 0.0, -1,
2876 	                     &dualvecsdiffer, NULL) );
2877 	         }
2878 	         else
2879 	         {
2880 	            SCIP_CALL( updateDualVector(scip, sepadata, dualvector, dualvector, nsoftcuts, 0, -1, -1, 0.0, NULL,
2881 	                     ngeneratedcurrroundcuts, TRUE, 0.0, 0.0, 0.0, 0.0, -1, &dualvecsdiffer, NULL) );
2882 	         }
2883 	      }
2884 	
2885 	      /* termination check */
2886 	      SCIP_CALL( checkMainLoopTermination(sepadata, cutoff, dualvecsdiffer, ngeneratedcurrroundcuts, nsoftcuts,
2887 	               nmaxgeneratedperroundcuts, ncurrroundlpiters, depth, &terminate) );
2888 	   }
2889 	
2890 	   /* set the maximal denominator in rational representation of gomory cut and the maximal scale factor to
2891 	    * scale resulting cut to integral values to avoid numerical instabilities
2892 	    */
2893 	   /**@todo find better but still stable gomory cut settings: look at dcmulti, gesa3, khb0525, misc06, p2756 */
2894 	   /* @note: above todo was copied from sepa_gomory.c. So, if gomory code is changed, same changes can be done here. */
2895 	   maxdepth = SCIPgetMaxDepth(scip);
2896 	   if( depth == 0 )
2897 	   {
2898 	      maxdnom = 1000;
2899 	      maxscale = 1000.0;
2900 	   }
2901 	   else if( depth <= maxdepth/4 )
2902 	   {
2903 	      maxdnom = 1000;
2904 	      maxscale = 1000.0;
2905 	   }
2906 	   else if( depth <= maxdepth/2 )
2907 	   {
2908 	      maxdnom = 100;
2909 	      maxscale = 100.0;
2910 	   }
2911 	   else
2912 	   {
2913 	      maxdnom = 10;
2914 	      maxscale = 10.0;
2915 	   }
2916 	
2917 	   /* filter cuts by calling cut selection algorithms and add cuts accordingly */
2918 	   /* @todo an idea is to filter cuts after every main iter */
2919 	   /* @todo we can remove !cutoff criterion by adding forcedcuts */
2920 	   if( !cutoff && (ngeneratedcurrroundcuts > 0) )
2921 	   {
2922 	      if( SCIPisGE(scip, sepadata->cutsfilterfactor, 1.0) )
2923 	      {
2924 	         nselectedcurrroundcuts = ngeneratedcurrroundcuts;
2925 	         SCIP_CALL( addCuts(scip, sepadata, generatedcurrroundcuts, nselectedcurrroundcuts, maxdnom, maxscale,
2926 	                  &naddedcurrroundcuts, &cutoff2) );
2927 	         cutoff = cutoff2;
2928 	      }
2929 	      else if( SCIPisPositive(scip, sepadata->cutsfilterfactor) )
2930 	      {
2931 	         nprocessedcuts = 0;
2932 	         for( int i = 0; i < sepadata->nmaxmainiters * (sepadata->nmaxsubgradientiters + 1); i++ )
2933 	         {
2934 	            if( ngeneratedcutsperiter[i] != 0 )
2935 	            {
2936 	               for( int j = 0; j < ngeneratedcutsperiter[i]; j++ )
2937 	                  cutindsperm[j] = j + nprocessedcuts;
2938 	
2939 	               /* sort cut efficacies by fractionality */
2940 	               SCIPsortDownRealInt(&generatedcutefficacies[nprocessedcuts], cutindsperm, ngeneratedcutsperiter[i]);
2941 	
2942 	               nselectedcuts = (int)SCIPceil(scip, sepadata->cutsfilterfactor * ngeneratedcutsperiter[i]);
2943 	
2944 	               SCIP_CALL( addCuts(scip, sepadata, &generatedcurrroundcuts[nprocessedcuts], nselectedcuts, maxdnom,
2945 	                        maxscale, &naddedcurrroundcuts, &cutoff2) );
2946 	               cutoff = cutoff2;
2947 	
2948 	               nprocessedcuts += ngeneratedcutsperiter[i];
2949 	            }
2950 	         }
2951 	      }
2952 	   }
2953 	   else if( ngeneratedcurrroundcuts > 0 )
2954 	   {
2955 	      nselectedcurrroundcuts = ngeneratedcurrroundcuts;
2956 	      SCIP_CALL( addCuts(scip, sepadata, generatedcurrroundcuts, nselectedcurrroundcuts, maxdnom, maxscale,
2957 	               &naddedcurrroundcuts, &cutoff2) );
2958 	   }
2959 	
2960 	   /* add an aggregated cut based on best Lagrangian multipliers */
2961 	   if( sepadata->aggregatecuts && (ngeneratedcurrroundcuts > 0) && (bestdualvectorlen > 0) )
2962 	   {
2963 	      assert(bestdualvectorlen <= ngeneratedcurrroundcuts);
2964 	      SCIP_CALL( aggregateGeneratedCuts(scip, sepa, sepadata, generatedcurrroundcuts, bestdualvector, bestdualvectorlen,
2965 	               aggregatedcurrroundcuts, &naggregatedcurrroundcuts, &cutoff2) );
2966 	      cutoff = (!cutoff ? cutoff2 : cutoff);
2967 	      if( naggregatedcurrroundcuts > 0 )
2968 	      {
2969 	         SCIP_CALL( addCuts(scip, sepadata, aggregatedcurrroundcuts, naggregatedcurrroundcuts, maxdnom, maxscale,
2970 	                  &naddedcurrroundcuts, &cutoff2) );
2971 	         cutoff = (!cutoff ? cutoff2 : cutoff);
2972 	      }
2973 	   }
2974 	
2975 	   if( cutoff )
2976 	   {
2977 	      *result = SCIP_CUTOFF;
2978 	   }
2979 	   else if( naddedcurrroundcuts > 0 )
2980 	   {
2981 	      *result = SCIP_SEPARATED;
2982 	   }
2983 	   else
2984 	   {
2985 	      *result = SCIP_DIDNOTFIND;
2986 	   }
2987 	
2988 	   /* delete the LP representing the Lagrangian dual problem with fixed Lagrangian multipliers */
2989 	   SCIP_CALL( deleteLPWithSoftCuts(scip, sepadata) );
2990 	
2991 	   /* release the rows in aggregatedcurrroundcuts */
2992 	   for( int i = 0; i < naggregatedcurrroundcuts; ++i )
2993 	   {
2994 	      SCIP_CALL( SCIPreleaseRow(scip, &(aggregatedcurrroundcuts[i])) );
2995 	   }
2996 	   /* release the rows in generatedcurrroundcuts */
2997 	   for( int i = 0; i < ngeneratedcurrroundcuts; ++i )
2998 	   {
2999 	      SCIP_CALL( SCIPreleaseRow(scip, &(generatedcurrroundcuts[i])) );
3000 	   }
3001 	
3002 	   for( int i = 0; i < sepadata->nmaxmainiters * (sepadata->nmaxsubgradientiters + 1); i++ )
3003 	   {
3004 	      ngeneratedcutsperiter[i] = 0;
3005 	   }
3006 	   for( int i = 0; i < nrows; i++ )
3007 	   {
3008 	      dualsol[i] = 0.0;
3009 	   }
3010 	   for( int i = 0; i < nmaxgeneratedperroundcuts; i++ )
3011 	   {
3012 	      dualsol[nrows + i] = 0.0;
3013 	      dualvector[i] = 0.0;
3014 	      bestdualvector[i] = 0.0;
3015 	   }
3016 	   /* free memory */
3017 	   SCIPfreeBufferArray(scip, &origobjcoefs);
3018 	   SCIPfreeBufferArray(scip, &hardcutslpsolvals);
3019 	   SCIP_CALL( SCIPfreeSol(scip, &hardcutslpsol) );
3020 	   SCIPfreeBufferArray(scip, &softcutslpsolvals);
3021 	   SCIP_CALL( SCIPfreeSol(scip, &softcutslpsol) );
3022 	   SCIPfreeCleanBufferArray(scip, &ngeneratedcutsperiter);
3023 	   SCIPfreeBufferArray(scip, &cutindsperm);
3024 	   SCIPfreeBufferArray(scip, &generatedcutefficacies);
3025 	   SCIPfreeBufferArray(scip, &aggregatedcurrroundcuts);
3026 	   SCIPfreeBufferArray(scip, &generatedcurrroundcuts);
3027 	   SCIPfreeCleanBufferArray(scip, &dualvector);
3028 	   SCIPfreeCleanBufferArray(scip, &bestdualvector);
3029 	   SCIPfreeCleanBufferArray(scip, &dualsol);
3030 	
3031 	   return SCIP_OKAY;
3032 	}
3033 	
3034 	/** creates separator data */
3035 	static
3036 	SCIP_RETCODE sepadataCreate(
3037 	   SCIP*                 scip,               /**< SCIP data structure */
3038 	   SCIP_SEPADATA**       sepadata            /**< separator data structure */
3039 	   )
3040 	{
3041 	   assert(scip != NULL);
3042 	   assert(sepadata != NULL);
3043 	
3044 	   SCIP_CALL( SCIPallocBlockMemory(scip, sepadata) );
3045 	   BMSclearMemory(*sepadata);
3046 	
3047 	   return SCIP_OKAY;
3048 	}
3049 	
3050 	
3051 	/*
3052 	 * Callback methods of separator
3053 	 */
3054 	
3055 	/** copy method for separator plugins (called when SCIP copies plugins) */
3056 	static
3057 	SCIP_DECL_SEPACOPY(sepaCopyLagromory)
3058 	{  /*lint --e{715}*/
3059 	   assert(scip != NULL);
3060 	   assert(sepa != NULL);
3061 	   assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
3062 	
3063 	   /* call inclusion method of constraint handler */
3064 	   SCIP_CALL( SCIPincludeSepaLagromory(scip) );
3065 	
3066 	   return SCIP_OKAY;
3067 	}
3068 	
3069 	
3070 	/** destructor of separator to free user data (called when SCIP is exiting) */
3071 	static
3072 	SCIP_DECL_SEPAFREE(sepaFreeLagromory)
3073 	{
3074 	   SCIP_SEPADATA* sepadata;
3075 	
3076 	   sepadata = SCIPsepaGetData(sepa);
3077 	   assert(sepadata != NULL);
3078 	
3079 	   SCIP_CALL( sepadataFree(scip, &sepadata) );
3080 	   SCIPsepaSetData(sepa, NULL);
3081 	
3082 	   return SCIP_OKAY;
3083 	}
3084 	
3085 	
3086 	/** initialization method of separator (called after problem was transformed) */
3087 	static
3088 	SCIP_DECL_SEPAINIT(sepaInitLagromory)
3089 	{
3090 	   SCIP_SEPADATA* sepadata;
3091 	
3092 	   sepadata = SCIPsepaGetData(sepa);
3093 	   assert(scip != NULL);
3094 	   assert(sepadata != NULL);
3095 	
3096 	   /* create and initialize random number generator */
3097 	   SCIP_CALL( SCIPcreateRandom(scip, &(sepadata->randnumgen), RANDSEED, TRUE) );
3098 	
3099 	   /* find trysol heuristic */
3100 	   if ( sepadata->heurtrysol == NULL )
3101 	   {
3102 	      sepadata->heurtrysol = SCIPfindHeur(scip, "trysol");
3103 	   }
3104 	
3105 	   return SCIP_OKAY;
3106 	}
3107 	
3108 	/** deinitialization method of separator (called before transformed problem is freed) */
3109 	static
3110 	SCIP_DECL_SEPAEXIT(sepaExitLagromory)
3111 	{  /*lint --e{715}*/
3112 	   SCIP_SEPADATA* sepadata;
3113 	
3114 	   sepadata = SCIPsepaGetData(sepa);
3115 	   assert(sepadata != NULL);
3116 	
3117 	   SCIPfreeRandom(scip, &(sepadata->randnumgen));
3118 	
3119 	   return SCIP_OKAY;
3120 	}
3121 	
3122 	/** LP solution separation method of separator */
3123 	static
3124 	SCIP_DECL_SEPAEXECLP(sepaExeclpLagromory)
3125 	{
3126 	   /*lint --e{715}*/
3127 	
3128 	   SCIP_SEPADATA* sepadata;
3129 	   SCIP_SOL* bestsol;
3130 	   SCIP_COL** cols;
3131 	   SCIP_COL* col;
3132 	   SCIP_VAR* var;
3133 	   SCIP_Real ubparam;
3134 	   SCIP_Real dualdegeneracyrate;
3135 	   SCIP_Real varconsratio;
3136 	   SCIP_Real threshold1;
3137 	   SCIP_Real threshold2;
3138 	   int nrows;
3139 	   int ncols;
3140 	   int ncalls;
3141 	   int runnum;
3142 	
3143 	   assert(sepa != NULL);
3144 	   assert(strcmp(SCIPsepaGetName(sepa), SEPA_NAME) == 0);
3145 	   sepadata = SCIPsepaGetData(sepa);
3146 	   assert(sepadata != NULL);
3147 	
3148 	   assert(result != NULL);
3149 	   *result = SCIP_DIDNOTRUN;
3150 	
3151 	   assert(scip != NULL);
3152 	
3153 	   ncalls = SCIPsepaGetNCallsAtNode(sepa);
3154 	   runnum = SCIPgetNRuns(scip);
3155 	   dualdegeneracyrate = 0.0;
3156 	   varconsratio = 0.0;
3157 	
3158 	   /* only generate Lagromory cuts if we are not close to terminating */
3159 	   if( SCIPisStopped(scip) )
3160 	      return SCIP_OKAY;
3161 	
3162 	   /* only call the separator starting sepadata->minrestart runs */
3163 	   if( (sepadata->minrestart >= 1) && (runnum < sepadata->minrestart + 1) )
3164 	      return SCIP_OKAY;
3165 	
3166 	   /* only call the separator a given number of times at each node */
3167 	   if( (depth == 0 && sepadata->maxroundsroot >= 0 && ncalls >= sepadata->maxroundsroot)
3168 	      || (depth > 0 && sepadata->maxrounds >= 0 && ncalls >= sepadata->maxrounds) )
3169 	      return SCIP_OKAY;
3170 	
3171 	   /* only generate cuts if an optimal LP solution is at hand */
3172 	   if( SCIPgetLPSolstat(scip) != SCIP_LPSOLSTAT_OPTIMAL )
3173 	      return SCIP_OKAY;
3174 	
3175 	   /* only generate cuts if the LP solution is basic */
3176 	   if( ! SCIPisLPSolBasic(scip) )
3177 	      return SCIP_OKAY;
3178 	
3179 	   /* get LP data */
3180 	   SCIP_CALL( SCIPgetLPColsData(scip, &cols, &ncols) );
3181 	   assert(cols != NULL);
3182 	
3183 	   nrows = SCIPgetNLPRows(scip);
3184 	
3185 	   /* return if LP has no columns or no rows */
3186 	   if( ncols == 0 || nrows == 0 )
3187 	      return SCIP_OKAY;
3188 	
3189 	   /* return if dual degeneracy metrics are below threshold values */
3190 	   threshold1 = sepadata->dualdegeneracyratethreshold;
3191 	   threshold2 = sepadata->varconsratiothreshold;
3192 	   SCIP_CALL( SCIPgetLPDualDegeneracy(scip, &dualdegeneracyrate, &varconsratio) );
3193 	   if( (dualdegeneracyrate < threshold1) && (varconsratio < threshold2) )
3194 	      return SCIP_OKAY;
3195 	
3196 	   bestsol = SCIPgetBestSol(scip);
3197 	   ubparam = 0.0;
3198 	
3199 	   /* if a MIP primal solution exists, use it to estimate the optimal value of the Lagrangian dual problem */
3200 	   if( bestsol != NULL )
3201 	   {
3202 	      for( int i = 0; i < ncols; ++i )
3203 	      {
3204 	         col = cols[i];
3205 	         assert(col != NULL);
3206 	
3207 	         var = SCIPcolGetVar(col);
3208 	
3209 	         ubparam += SCIPgetSolVal(scip, bestsol, var) * SCIPcolGetObj(col);
3210 	      }
3211 	
3212 	      ubparam += SCIPgetTransObjoffset(scip);
3213 	   }
3214 	   /* if a MIP primal solution does not exist, then use the node relaxation's LP solutin to estimate the optimal value
3215 	    * of the Lagrangian dual problem
3216 	    */
3217 	   else
3218 	   {
3219 	      for( int i = 0; i < ncols; ++i )
3220 	      {
3221 	         col = cols[i];
3222 	         assert(col != NULL);
3223 	
3224 	         ubparam += SCIPcolGetPrimsol(col) * SCIPcolGetObj(col);
3225 	      }
3226 	
3227 	      ubparam += SCIPgetTransObjoffset(scip);
3228 	      ubparam *= SCIPisPositive(scip, ubparam) ? sepadata->ubparamposfactor : sepadata->ubparamnegfactor;
3229 	   }
3230 	
3231 	   /* the main separation function of the separator */
3232 	   SCIP_CALL( separateCuts(scip, sepa, sepadata, ubparam, depth, (allowlocal && sepadata->allowlocal), result) );
3233 	
3234 	   return SCIP_OKAY;
3235 	}
3236 	
3237 	
3238 	/*
3239 	 * separator specific interface methods
3240 	 */
3241 	
3242 	/** creates the Lagromory separator and includes it in SCIP */
3243 	SCIP_RETCODE SCIPincludeSepaLagromory(
3244 	   SCIP*                 scip                /**< SCIP data structure */
3245 	   )
3246 	{
3247 	   SCIP_SEPADATA* sepadata;
3248 	   SCIP_SEPA* sepa;
3249 	
3250 	   /* create separator data */
3251 	   SCIP_CALL( sepadataCreate(scip, &sepadata) );
3252 	
3253 	   sepadata->heurtrysol = NULL;
3254 	
3255 	   /* include separator */
3256 	   SCIP_CALL( SCIPincludeSepaBasic(scip, &sepa, SEPA_NAME, SEPA_DESC, SEPA_PRIORITY, SEPA_FREQ, SEPA_MAXBOUNDDIST,
3257 	            SEPA_USESSUBSCIP, SEPA_DELAY, sepaExeclpLagromory, NULL, sepadata) );
3258 	
3259 	   assert(sepa != NULL);
3260 	
3261 	   /* set non fundamental callbacks via setter functions */
3262 	   SCIP_CALL( SCIPsetSepaCopy(scip, sepa, sepaCopyLagromory) );
3263 	   SCIP_CALL( SCIPsetSepaFree(scip, sepa, sepaFreeLagromory) );
3264 	   SCIP_CALL( SCIPsetSepaInit(scip, sepa, sepaInitLagromory) );
3265 	   SCIP_CALL( SCIPsetSepaExit(scip, sepa, sepaExitLagromory) );
3266 	
3267 	   /* add separator parameters */
3268 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/away", "minimal integrality violation of a basis "
3269 	            "variable to try separation", &sepadata->away, FALSE, DEFAULT_AWAY, 0.0, 1.0, NULL, NULL) );
3270 	
3271 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/rootlpiterlimitfactor", "factor w.r.t. root node LP "
3272 	            "iterations for maximal separating LP iterations in the root node (negative for no limit)",
3273 	            &sepadata->rootlpiterlimitfactor, TRUE, DEFAULT_ROOTLPITERLIMITFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
3274 	
3275 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/totallpiterlimitfactor", "factor w.r.t. root node LP "
3276 	            "iterations for maximal separating LP iterations in the tree (negative for no limit)",
3277 	            &sepadata->totallpiterlimitfactor, TRUE, DEFAULT_TOTALLPITERLIMITFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
3278 	
3279 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/perroundlpiterlimitfactor", "factor w.r.t. root node LP "
3280 	            "iterations for maximal separating LP iterations per separation round (negative for no limit)",
3281 	            &sepadata->perroundlpiterlimitfactor, TRUE, DEFAULT_PERROUNDLPITERLIMITFACTOR, -1.0, SCIP_REAL_MAX, NULL,
3282 	            NULL) );
3283 	
3284 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/perroundcutsfactorroot", "factor w.r.t. number of integer "
3285 	            "columns for number of cuts separated per separation round in root node", &sepadata->perroundcutsfactorroot,
3286 	            TRUE, DEFAULT_PERROUNDCUTSFACTORROOT, 0.0, SCIP_REAL_MAX, NULL, NULL) );
3287 	
3288 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/perroundcutsfactor", "factor w.r.t. number of integer "
3289 	            "columns for number of cuts separated per separation round at a non-root node",
3290 	            &sepadata->perroundcutsfactor, TRUE, DEFAULT_PERROUNDCUTSFACTOR, 0.0, SCIP_REAL_MAX, NULL, NULL) );
3291 	
3292 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/totalcutsfactor", "factor w.r.t. number of integer "
3293 	            "columns for total number of cuts separated", &sepadata->totalcutsfactor, TRUE, DEFAULT_TOTALCUTSFACTOR, 0.0,
3294 	            SCIP_REAL_MAX, NULL, NULL) );
3295 	
3296 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muparaminit", "initial value of the mu parameter (factor "
3297 	      "for step length)", &sepadata->muparaminit, TRUE, DEFAULT_MUPARAMINIT, 0.0, 100.0, NULL, NULL) );
3298 	
3299 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muparamlb", "lower bound of the mu parameter (factor for "
3300 	      "step length)", &sepadata->muparamlb, TRUE, DEFAULT_MUPARAMLB, 0.0, 1.0, NULL, NULL) );
3301 	
3302 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muparamub", "upper bound of the mu parameter (factor for "
3303 	      "step length)", &sepadata->muparamub, TRUE, DEFAULT_MUPARAMUB, 1.0, 10.0, NULL, NULL) );
3304 	
3305 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/mubacktrackfactor", "factor of mu while backtracking the "
3306 	            "mu parameter (factor for step length)", &sepadata->mubacktrackfactor, TRUE, DEFAULT_MUBACKTRACKFACTOR,
3307 	            0.0, 1.0, NULL, NULL) );
3308 	
3309 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muslab1factor", "factor of mu parameter (factor for step "
3310 	      "length) for larger increment" , &sepadata->muslab1factor, TRUE, DEFAULT_MUSLAB1FACTOR, 0.0, SCIP_REAL_MAX, NULL,
3311 	            NULL));
3312 	
3313 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muslab2factor", "factor of mu parameter (factor for step "
3314 	      "length) for smaller increment", &sepadata->muslab2factor, TRUE, DEFAULT_MUSLAB2FACTOR, 0.0, SCIP_REAL_MAX, NULL,
3315 	            NULL) );
3316 	
3317 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/muslab3factor", "factor of mu parameter (factor for step "
3318 	      "length) for reduction", &sepadata->muslab3factor, TRUE, DEFAULT_MUSLAB3FACTOR, 0.0, SCIP_REAL_MAX, NULL, NULL) );
3319 	
3320 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/deltaslab1ub", "factor of delta deciding larger increment "
3321 	            "of mu parameter (factor for step length)", &sepadata->deltaslab1ub, TRUE, DEFAULT_DELTASLAB1UB, 0.0, 1.0,
3322 	            NULL, NULL) );
3323 	
3324 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/deltaslab2ub", "factor of delta deciding smaller "
3325 	            "increment of mu parameter (factor for step length)", &sepadata->deltaslab2ub, TRUE, DEFAULT_DELTASLAB2UB,
3326 	            0.0, 1.0, NULL, NULL) );
3327 	
3328 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/ubparamposfactor", "factor for positive upper bound used "
3329 	            "as an estimate for the optimal Lagrangian dual value", &sepadata->ubparamposfactor, TRUE, DEFAULT_UBPARAMPOSFACTOR,
3330 	            1.0, 100.0, NULL, NULL) );
3331 	
3332 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/ubparamnegfactor", "factor for negative upper bound used "
3333 	            "as an estimate for the optimal Lagrangian dual value", &sepadata->ubparamnegfactor, TRUE, DEFAULT_UBPARAMNEGFACTOR,
3334 	            0.0, 1.0, NULL, NULL) );
3335 	
3336 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/perrootlpiterfactor", "factor w.r.t. root node LP "
3337 	            "iterations for iteration limit of each separating LP (negative for no limit)",
3338 	            &sepadata->perrootlpiterfactor, TRUE, DEFAULT_PERROOTLPITERFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
3339 	
3340 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/perlpiterfactor", "factor w.r.t. non-root node LP "
3341 	            "iterations for iteration limit of each separating LP (negative for no limit)", &sepadata->perlpiterfactor, TRUE,
3342 	            DEFAULT_PERLPITERFACTOR, -1.0, SCIP_REAL_MAX, NULL, NULL) );
3343 	
3344 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/cutsfilterfactor", "fraction of generated cuts per "
3345 	            "explored basis to accept from separator", &sepadata->cutsfilterfactor, TRUE, DEFAULT_CUTSFILTERFACTOR, 0.0,
3346 	            1.0, NULL, NULL) );
3347 	
3348 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/radiusinit", "initial radius of the ball used in "
3349 	            "stabilization of Lagrangian multipliers", &sepadata->radiusinit, TRUE, DEFAULT_RADIUSINIT, 0.0, 1.0, NULL,
3350 	            NULL) );
3351 	
3352 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/radiusmax", "maximum radius of the ball used in "
3353 	            "stabilization of Lagrangian multipliers", &sepadata->radiusmax, TRUE, DEFAULT_RADIUSMAX, 0.0, SCIP_REAL_MAX,
3354 	            NULL, NULL) );
3355 	
3356 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/radiusmin", "minimum radius of the ball used in "
3357 	            "stabilization of Lagrangian multipliers", &sepadata->radiusmin, TRUE, DEFAULT_RADIUSMIN, 0.0, SCIP_REAL_MAX,
3358 	            NULL, NULL) );
3359 	
3360 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/constant", "a constant for stablity center based "
3361 	            "stabilization of Lagrangian multipliers", &sepadata->constant, TRUE, DEFAULT_CONST, 2.0, SCIP_REAL_MAX,
3362 	            NULL, NULL) );
3363 	
3364 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/radiusupdateweight", "multiplier to evaluate cut "
3365 	            "violation score used for updating ball radius", &sepadata->radiusupdateweight, TRUE,
3366 	            DEFAULT_RADIUSUPDATEWEIGHT, 0.0, 1.0, NULL, NULL) );
3367 	
3368 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/dualdegeneracyratethreshold", "minimum dual degeneracy "
3369 	            "rate for separator execution", &sepadata->dualdegeneracyratethreshold, FALSE,
3370 	            DEFAULT_DUALDEGENERACYRATETHRESHOLD, 0.0, 1.0, NULL, NULL) );
3371 	
3372 	   SCIP_CALL( SCIPaddRealParam(scip, "separating/" SEPA_NAME "/varconsratiothreshold", "minimum variable-constraint "
3373 	            "ratio on optimal face for separator execution", &sepadata->varconsratiothreshold, FALSE,
3374 	            DEFAULT_VARCONSRATIOTHRESHOLD, 1.0, SCIP_REAL_MAX, NULL, NULL) );
3375 	
3376 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/muparamconst", "is the mu parameter (factor for step "
3377 	      "length) constant?" , &sepadata->muparamconst, TRUE, DEFAULT_MUPARAMCONST, NULL, NULL) );
3378 	
3379 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/separaterows", "separate rows with integral slack?",
3380 	            &sepadata->separaterows, TRUE, DEFAULT_SEPARATEROWS, NULL, NULL) );
3381 	
3382 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/sortcutoffsol", "sort fractional integer columns"
3383 	            "based on fractionality?" , &sepadata->sortcutoffsol, TRUE, DEFAULT_SORTCUTOFFSOL, NULL, NULL) );
3384 	
3385 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/sidetypebasis", "choose side types of row (lhs/rhs) "
3386 	            "based on basis information?", &sepadata->sidetypebasis, TRUE, DEFAULT_SIDETYPEBASIS, NULL, NULL) );
3387 	
3388 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/dynamiccuts", "should generated cuts be removed from "
3389 	            "LP if they are no longer tight?", &sepadata->dynamiccuts, FALSE, DEFAULT_DYNAMICCUTS, NULL, NULL) );
3390 	
3391 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/makeintegral", "try to scale all cuts to integral "
3392 	            "coefficients?", &sepadata->makeintegral, TRUE, DEFAULT_MAKEINTEGRAL, NULL, NULL) );
3393 	
3394 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/forcecuts", "force cuts to be added to the LP?",
3395 	            &sepadata->forcecuts, TRUE, DEFAULT_FORCECUTS, NULL, NULL) );
3396 	
3397 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/delayedcuts", "should cuts be added to the delayed cut "
3398 	            "pool", &sepadata->delayedcuts, TRUE, DEFAULT_DELAYEDCUTS, NULL, NULL) );
3399 	
3400 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/allowlocal", "should locally valid cuts be generated?",
3401 	            &sepadata->allowlocal, TRUE, DEFAULT_ALLOWLOCAL, NULL, NULL) );
3402 	
3403 	   SCIP_CALL( SCIPaddBoolParam(scip, "separating/" SEPA_NAME "/aggregatecuts", "aggregate all generated cuts using the "
3404 	            "Lagrangian multipliers?", &sepadata->aggregatecuts, TRUE, DEFAULT_AGGREGATECUTS, NULL, NULL) );
3405 	
3406 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/maxrounds", "maximal number of separation rounds per node "
3407 	            "(-1: unlimited)", &sepadata->maxrounds, FALSE, DEFAULT_MAXROUNDS, -1, INT_MAX, NULL, NULL) );
3408 	
3409 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/maxroundsroot", "maximal number of separation rounds in "
3410 	            "the root node (-1: unlimited)", &sepadata->maxroundsroot, FALSE, DEFAULT_MAXROUNDSROOT, -1, INT_MAX, NULL,
3411 	            NULL) );
3412 	
3413 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/perroundnmaxlpiters", "maximal number of separating LP "
3414 	            "iterations per separation round (-1: unlimited)", &sepadata->perroundnmaxlpiters, FALSE,
3415 	            DEFAULT_PERROUNDMAXLPITERS, -1, INT_MAX, NULL, NULL) );
3416 	
3417 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxcutsperlp", "maximal number of cuts separated per "
3418 	            "Lagromory LP in the non-root node", &sepadata->nmaxcutsperlp, FALSE, DEFAULT_PERLPMAXCUTS, 0, INT_MAX, NULL,
3419 	            NULL));
3420 	
3421 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxcutsperlproot", "maximal number of cuts separated per "
3422 	            "Lagromory LP in the root node", &sepadata->nmaxcutsperlproot, FALSE, DEFAULT_PERLPMAXCUTSROOT, 0, INT_MAX,
3423 	            NULL, NULL));
3424 	
3425 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxmainiters", "maximal number of main loop iterations of "
3426 	            "the relax-and-cut algorithm", &sepadata->nmaxmainiters, TRUE, DEFAULT_MAXMAINITERS, 0, INT_MAX, NULL, NULL)
3427 	         );
3428 	
3429 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxsubgradientiters", "maximal number of subgradient loop "
3430 	            "iterations of the relax-and-cut algorithm", &sepadata->nmaxsubgradientiters, TRUE,
3431 	            DEFAULT_MAXSUBGRADIENTITERS, 0, INT_MAX, NULL, NULL) );
3432 	
3433 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/cutgenfreq", "frequency of subgradient iterations for "
3434 	            "generating cuts", &sepadata->cutgenfreq, TRUE, DEFAULT_CUTGENFREQ, 0, INT_MAX, NULL, NULL) );
3435 	
3436 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/cutaddfreq", "frequency of subgradient iterations for "
3437 	            "adding cuts to objective function", &sepadata->cutaddfreq, TRUE, DEFAULT_CUTADDFREQ, 0, INT_MAX, NULL, NULL) );
3438 	
3439 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxlagrangianvalsforavg", "maximal number of iterations "
3440 	            "for rolling average of Lagrangian value", &sepadata->nmaxlagrangianvalsforavg, TRUE,
3441 	            DEFAULT_MAXLAGRANGIANVALSFORAVG, 0, INT_MAX, NULL, NULL) );
3442 	
3443 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/nmaxconsecitersformuupdate", "consecutive number of "
3444 	            "iterations used to determine if mu needs to be backtracked", &sepadata->nmaxconsecitersformuupdate, TRUE,
3445 	            DEFAULT_MAXCONSECITERSFORMUUPDATE, 0, INT_MAX, NULL, NULL) );
3446 	
3447 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/projectiontype", "the ball into which the Lagrangian multipliers "
3448 	            "are projected for stabilization (0: no projection, 1: L1-norm ball projection, 2: L2-norm ball projection, 3: "
3449 	            "L_inf-norm ball projection)", &sepadata->projectiontype, TRUE, DEFAULT_PROJECTIONTYPE, 0, 3, NULL, NULL));
3450 	
3451 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/stabilitycentertype", "type of stability center for "
3452 	            "taking weighted average of Lagrangian multipliers for stabilization (0: no weighted stabilization, 1: best "
3453 	            "Lagrangian multipliers)", &sepadata->stabilitycentertype, TRUE, DEFAULT_STABILITYCENTERTYPE, 0, 1, NULL,
3454 	            NULL) );
3455 	
3456 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/optimalfacepriority", "priority of the optimal face for "
3457 	            "separator execution (0: low priority, 1: medium priority, 2: high priority)",
3458 	            &sepadata->optimalfacepriority, TRUE, DEFAULT_OPTIMALFACEPRIORITY, 0, 2, NULL, NULL) );
3459 	
3460 	   SCIP_CALL( SCIPaddIntParam(scip, "separating/" SEPA_NAME "/minrestart", "minimum restart round for separator "
3461 	            "execution (0: from beginning of the instance solving, >= n with n >= 1: from restart round n)",
3462 	            &sepadata->minrestart, TRUE, DEFAULT_MINRESTART, 0, INT_MAX, NULL, NULL) );
3463 	
3464 	  return SCIP_OKAY;
3465 	}
3466