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