Scippy

SCIP

Solving Constraint Integer Programs

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