1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/glpios03.c Mon Dec 06 13:09:21 2010 +0100
1.3 @@ -0,0 +1,1155 @@
1.4 +/* glpios03.c (branch-and-cut driver) */
1.5 +
1.6 +/***********************************************************************
1.7 +* This code is part of GLPK (GNU Linear Programming Kit).
1.8 +*
1.9 +* Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
1.10 +* 2009, 2010 Andrew Makhorin, Department for Applied Informatics,
1.11 +* Moscow Aviation Institute, Moscow, Russia. All rights reserved.
1.12 +* E-mail: <mao@gnu.org>.
1.13 +*
1.14 +* GLPK is free software: you can redistribute it and/or modify it
1.15 +* under the terms of the GNU General Public License as published by
1.16 +* the Free Software Foundation, either version 3 of the License, or
1.17 +* (at your option) any later version.
1.18 +*
1.19 +* GLPK is distributed in the hope that it will be useful, but WITHOUT
1.20 +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.21 +* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
1.22 +* License for more details.
1.23 +*
1.24 +* You should have received a copy of the GNU General Public License
1.25 +* along with GLPK. If not, see <http://www.gnu.org/licenses/>.
1.26 +***********************************************************************/
1.27 +
1.28 +#include "glpios.h"
1.29 +
1.30 +/***********************************************************************
1.31 +* show_progress - display current progress of the search
1.32 +*
1.33 +* This routine displays some information about current progress of the
1.34 +* search.
1.35 +*
1.36 +* The information includes:
1.37 +*
1.38 +* the current number of iterations performed by the simplex solver;
1.39 +*
1.40 +* the objective value for the best known integer feasible solution,
1.41 +* which is upper (minimization) or lower (maximization) global bound
1.42 +* for optimal solution of the original mip problem;
1.43 +*
1.44 +* the best local bound for active nodes, which is lower (minimization)
1.45 +* or upper (maximization) global bound for optimal solution of the
1.46 +* original mip problem;
1.47 +*
1.48 +* the relative mip gap, in percents;
1.49 +*
1.50 +* the number of open (active) subproblems;
1.51 +*
1.52 +* the number of completely explored subproblems, i.e. whose nodes have
1.53 +* been removed from the tree. */
1.54 +
1.55 +static void show_progress(glp_tree *T, int bingo)
1.56 +{ int p;
1.57 + double temp;
1.58 + char best_mip[50], best_bound[50], *rho, rel_gap[50];
1.59 + /* format the best known integer feasible solution */
1.60 + if (T->mip->mip_stat == GLP_FEAS)
1.61 + sprintf(best_mip, "%17.9e", T->mip->mip_obj);
1.62 + else
1.63 + sprintf(best_mip, "%17s", "not found yet");
1.64 + /* determine reference number of an active subproblem whose local
1.65 + bound is best */
1.66 + p = ios_best_node(T);
1.67 + /* format the best bound */
1.68 + if (p == 0)
1.69 + sprintf(best_bound, "%17s", "tree is empty");
1.70 + else
1.71 + { temp = T->slot[p].node->bound;
1.72 + if (temp == -DBL_MAX)
1.73 + sprintf(best_bound, "%17s", "-inf");
1.74 + else if (temp == +DBL_MAX)
1.75 + sprintf(best_bound, "%17s", "+inf");
1.76 + else
1.77 + sprintf(best_bound, "%17.9e", temp);
1.78 + }
1.79 + /* choose the relation sign between global bounds */
1.80 + if (T->mip->dir == GLP_MIN)
1.81 + rho = ">=";
1.82 + else if (T->mip->dir == GLP_MAX)
1.83 + rho = "<=";
1.84 + else
1.85 + xassert(T != T);
1.86 + /* format the relative mip gap */
1.87 + temp = ios_relative_gap(T);
1.88 + if (temp == 0.0)
1.89 + sprintf(rel_gap, " 0.0%%");
1.90 + else if (temp < 0.001)
1.91 + sprintf(rel_gap, "< 0.1%%");
1.92 + else if (temp <= 9.999)
1.93 + sprintf(rel_gap, "%5.1f%%", 100.0 * temp);
1.94 + else
1.95 + sprintf(rel_gap, "%6s", "");
1.96 + /* display progress of the search */
1.97 + xprintf("+%6d: %s %s %s %s %s (%d; %d)\n",
1.98 + T->mip->it_cnt, bingo ? ">>>>>" : "mip =", best_mip, rho,
1.99 + best_bound, rel_gap, T->a_cnt, T->t_cnt - T->n_cnt);
1.100 + T->tm_lag = xtime();
1.101 + return;
1.102 +}
1.103 +
1.104 +/***********************************************************************
1.105 +* is_branch_hopeful - check if specified branch is hopeful
1.106 +*
1.107 +* This routine checks if the specified subproblem can have an integer
1.108 +* optimal solution which is better than the best known one.
1.109 +*
1.110 +* The check is based on comparison of the local objective bound stored
1.111 +* in the subproblem descriptor and the incumbent objective value which
1.112 +* is the global objective bound.
1.113 +*
1.114 +* If there is a chance that the specified subproblem can have a better
1.115 +* integer optimal solution, the routine returns non-zero. Otherwise, if
1.116 +* the corresponding branch can pruned, zero is returned. */
1.117 +
1.118 +static int is_branch_hopeful(glp_tree *T, int p)
1.119 +{ xassert(1 <= p && p <= T->nslots);
1.120 + xassert(T->slot[p].node != NULL);
1.121 + return ios_is_hopeful(T, T->slot[p].node->bound);
1.122 +}
1.123 +
1.124 +/***********************************************************************
1.125 +* check_integrality - check integrality of basic solution
1.126 +*
1.127 +* This routine checks if the basic solution of LP relaxation of the
1.128 +* current subproblem satisfies to integrality conditions, i.e. that all
1.129 +* variables of integer kind have integral primal values. (The solution
1.130 +* is assumed to be optimal.)
1.131 +*
1.132 +* For each variable of integer kind the routine computes the following
1.133 +* quantity:
1.134 +*
1.135 +* ii(x[j]) = min(x[j] - floor(x[j]), ceil(x[j]) - x[j]), (1)
1.136 +*
1.137 +* which is a measure of the integer infeasibility (non-integrality) of
1.138 +* x[j] (for example, ii(2.1) = 0.1, ii(3.7) = 0.3, ii(5.0) = 0). It is
1.139 +* understood that 0 <= ii(x[j]) <= 0.5, and variable x[j] is integer
1.140 +* feasible if ii(x[j]) = 0. However, due to floating-point arithmetic
1.141 +* the routine checks less restrictive condition:
1.142 +*
1.143 +* ii(x[j]) <= tol_int, (2)
1.144 +*
1.145 +* where tol_int is a given tolerance (small positive number) and marks
1.146 +* each variable which does not satisfy to (2) as integer infeasible by
1.147 +* setting its fractionality flag.
1.148 +*
1.149 +* In order to characterize integer infeasibility of the basic solution
1.150 +* in the whole the routine computes two parameters: ii_cnt, which is
1.151 +* the number of variables with the fractionality flag set, and ii_sum,
1.152 +* which is the sum of integer infeasibilities (1). */
1.153 +
1.154 +static void check_integrality(glp_tree *T)
1.155 +{ glp_prob *mip = T->mip;
1.156 + int j, type, ii_cnt = 0;
1.157 + double lb, ub, x, temp1, temp2, ii_sum = 0.0;
1.158 + /* walk through the set of columns (structural variables) */
1.159 + for (j = 1; j <= mip->n; j++)
1.160 + { GLPCOL *col = mip->col[j];
1.161 + T->non_int[j] = 0;
1.162 + /* if the column is not integer, skip it */
1.163 + if (col->kind != GLP_IV) continue;
1.164 + /* if the column is non-basic, it is integer feasible */
1.165 + if (col->stat != GLP_BS) continue;
1.166 + /* obtain the type and bounds of the column */
1.167 + type = col->type, lb = col->lb, ub = col->ub;
1.168 + /* obtain value of the column in optimal basic solution */
1.169 + x = col->prim;
1.170 + /* if the column's primal value is close to the lower bound,
1.171 + the column is integer feasible within given tolerance */
1.172 + if (type == GLP_LO || type == GLP_DB || type == GLP_FX)
1.173 + { temp1 = lb - T->parm->tol_int;
1.174 + temp2 = lb + T->parm->tol_int;
1.175 + if (temp1 <= x && x <= temp2) continue;
1.176 +#if 0
1.177 + /* the lower bound must not be violated */
1.178 + xassert(x >= lb);
1.179 +#else
1.180 + if (x < lb) continue;
1.181 +#endif
1.182 + }
1.183 + /* if the column's primal value is close to the upper bound,
1.184 + the column is integer feasible within given tolerance */
1.185 + if (type == GLP_UP || type == GLP_DB || type == GLP_FX)
1.186 + { temp1 = ub - T->parm->tol_int;
1.187 + temp2 = ub + T->parm->tol_int;
1.188 + if (temp1 <= x && x <= temp2) continue;
1.189 +#if 0
1.190 + /* the upper bound must not be violated */
1.191 + xassert(x <= ub);
1.192 +#else
1.193 + if (x > ub) continue;
1.194 +#endif
1.195 + }
1.196 + /* if the column's primal value is close to nearest integer,
1.197 + the column is integer feasible within given tolerance */
1.198 + temp1 = floor(x + 0.5) - T->parm->tol_int;
1.199 + temp2 = floor(x + 0.5) + T->parm->tol_int;
1.200 + if (temp1 <= x && x <= temp2) continue;
1.201 + /* otherwise the column is integer infeasible */
1.202 + T->non_int[j] = 1;
1.203 + /* increase the number of fractional-valued columns */
1.204 + ii_cnt++;
1.205 + /* compute the sum of integer infeasibilities */
1.206 + temp1 = x - floor(x);
1.207 + temp2 = ceil(x) - x;
1.208 + xassert(temp1 > 0.0 && temp2 > 0.0);
1.209 + ii_sum += (temp1 <= temp2 ? temp1 : temp2);
1.210 + }
1.211 + /* store ii_cnt and ii_sum to the current problem descriptor */
1.212 + xassert(T->curr != NULL);
1.213 + T->curr->ii_cnt = ii_cnt;
1.214 + T->curr->ii_sum = ii_sum;
1.215 + /* and also display these parameters */
1.216 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.217 + { if (ii_cnt == 0)
1.218 + xprintf("There are no fractional columns\n");
1.219 + else if (ii_cnt == 1)
1.220 + xprintf("There is one fractional column, integer infeasibil"
1.221 + "ity is %.3e\n", ii_sum);
1.222 + else
1.223 + xprintf("There are %d fractional columns, integer infeasibi"
1.224 + "lity is %.3e\n", ii_cnt, ii_sum);
1.225 + }
1.226 + return;
1.227 +}
1.228 +
1.229 +/***********************************************************************
1.230 +* record_solution - record better integer feasible solution
1.231 +*
1.232 +* This routine records optimal basic solution of LP relaxation of the
1.233 +* current subproblem, which being integer feasible is better than the
1.234 +* best known integer feasible solution. */
1.235 +
1.236 +static void record_solution(glp_tree *T)
1.237 +{ glp_prob *mip = T->mip;
1.238 + int i, j;
1.239 + mip->mip_stat = GLP_FEAS;
1.240 + mip->mip_obj = mip->obj_val;
1.241 + for (i = 1; i <= mip->m; i++)
1.242 + { GLPROW *row = mip->row[i];
1.243 + row->mipx = row->prim;
1.244 + }
1.245 + for (j = 1; j <= mip->n; j++)
1.246 + { GLPCOL *col = mip->col[j];
1.247 + if (col->kind == GLP_CV)
1.248 + col->mipx = col->prim;
1.249 + else if (col->kind == GLP_IV)
1.250 + { /* value of the integer column must be integral */
1.251 + col->mipx = floor(col->prim + 0.5);
1.252 + }
1.253 + else
1.254 + xassert(col != col);
1.255 + }
1.256 + T->sol_cnt++;
1.257 + return;
1.258 +}
1.259 +
1.260 +/***********************************************************************
1.261 +* fix_by_red_cost - fix non-basic integer columns by reduced costs
1.262 +*
1.263 +* This routine fixes some non-basic integer columns if their reduced
1.264 +* costs indicate that increasing (decreasing) the column at least by
1.265 +* one involves the objective value becoming worse than the incumbent
1.266 +* objective value. */
1.267 +
1.268 +static void fix_by_red_cost(glp_tree *T)
1.269 +{ glp_prob *mip = T->mip;
1.270 + int j, stat, fixed = 0;
1.271 + double obj, lb, ub, dj;
1.272 + /* the global bound must exist */
1.273 + xassert(T->mip->mip_stat == GLP_FEAS);
1.274 + /* basic solution of LP relaxation must be optimal */
1.275 + xassert(mip->pbs_stat == GLP_FEAS && mip->dbs_stat == GLP_FEAS);
1.276 + /* determine the objective function value */
1.277 + obj = mip->obj_val;
1.278 + /* walk through the column list */
1.279 + for (j = 1; j <= mip->n; j++)
1.280 + { GLPCOL *col = mip->col[j];
1.281 + /* if the column is not integer, skip it */
1.282 + if (col->kind != GLP_IV) continue;
1.283 + /* obtain bounds of j-th column */
1.284 + lb = col->lb, ub = col->ub;
1.285 + /* and determine its status and reduced cost */
1.286 + stat = col->stat, dj = col->dual;
1.287 + /* analyze the reduced cost */
1.288 + switch (mip->dir)
1.289 + { case GLP_MIN:
1.290 + /* minimization */
1.291 + if (stat == GLP_NL)
1.292 + { /* j-th column is non-basic on its lower bound */
1.293 + if (dj < 0.0) dj = 0.0;
1.294 + if (obj + dj >= mip->mip_obj)
1.295 + glp_set_col_bnds(mip, j, GLP_FX, lb, lb), fixed++;
1.296 + }
1.297 + else if (stat == GLP_NU)
1.298 + { /* j-th column is non-basic on its upper bound */
1.299 + if (dj > 0.0) dj = 0.0;
1.300 + if (obj - dj >= mip->mip_obj)
1.301 + glp_set_col_bnds(mip, j, GLP_FX, ub, ub), fixed++;
1.302 + }
1.303 + break;
1.304 + case GLP_MAX:
1.305 + /* maximization */
1.306 + if (stat == GLP_NL)
1.307 + { /* j-th column is non-basic on its lower bound */
1.308 + if (dj > 0.0) dj = 0.0;
1.309 + if (obj + dj <= mip->mip_obj)
1.310 + glp_set_col_bnds(mip, j, GLP_FX, lb, lb), fixed++;
1.311 + }
1.312 + else if (stat == GLP_NU)
1.313 + { /* j-th column is non-basic on its upper bound */
1.314 + if (dj < 0.0) dj = 0.0;
1.315 + if (obj - dj <= mip->mip_obj)
1.316 + glp_set_col_bnds(mip, j, GLP_FX, ub, ub), fixed++;
1.317 + }
1.318 + break;
1.319 + default:
1.320 + xassert(T != T);
1.321 + }
1.322 + }
1.323 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.324 + { if (fixed == 0)
1.325 + /* nothing to say */;
1.326 + else if (fixed == 1)
1.327 + xprintf("One column has been fixed by reduced cost\n");
1.328 + else
1.329 + xprintf("%d columns have been fixed by reduced costs\n",
1.330 + fixed);
1.331 + }
1.332 + /* fixing non-basic columns on their current bounds does not
1.333 + change the basic solution */
1.334 + xassert(mip->pbs_stat == GLP_FEAS && mip->dbs_stat == GLP_FEAS);
1.335 + return;
1.336 +}
1.337 +
1.338 +/***********************************************************************
1.339 +* branch_on - perform branching on specified variable
1.340 +*
1.341 +* This routine performs branching on j-th column (structural variable)
1.342 +* of the current subproblem. The specified column must be of integer
1.343 +* kind and must have a fractional value in optimal basic solution of
1.344 +* LP relaxation of the current subproblem (i.e. only columns for which
1.345 +* the flag non_int[j] is set are valid candidates to branch on).
1.346 +*
1.347 +* Let x be j-th structural variable, and beta be its primal fractional
1.348 +* value in the current basic solution. Branching on j-th variable is
1.349 +* dividing the current subproblem into two new subproblems, which are
1.350 +* identical to the current subproblem with the following exception: in
1.351 +* the first subproblem that begins the down-branch x has a new upper
1.352 +* bound x <= floor(beta), and in the second subproblem that begins the
1.353 +* up-branch x has a new lower bound x >= ceil(beta).
1.354 +*
1.355 +* Depending on estimation of local bounds for down- and up-branches
1.356 +* this routine returns the following:
1.357 +*
1.358 +* 0 - both branches have been created;
1.359 +* 1 - one branch is hopeless and has been pruned, so now the current
1.360 +* subproblem is other branch;
1.361 +* 2 - both branches are hopeless and have been pruned; new subproblem
1.362 +* selection is needed to continue the search. */
1.363 +
1.364 +static int branch_on(glp_tree *T, int j, int next)
1.365 +{ glp_prob *mip = T->mip;
1.366 + IOSNPD *node;
1.367 + int m = mip->m;
1.368 + int n = mip->n;
1.369 + int type, dn_type, up_type, dn_bad, up_bad, p, ret, clone[1+2];
1.370 + double lb, ub, beta, new_ub, new_lb, dn_lp, up_lp, dn_bnd, up_bnd;
1.371 + /* determine bounds and value of x[j] in optimal solution to LP
1.372 + relaxation of the current subproblem */
1.373 + xassert(1 <= j && j <= n);
1.374 + type = mip->col[j]->type;
1.375 + lb = mip->col[j]->lb;
1.376 + ub = mip->col[j]->ub;
1.377 + beta = mip->col[j]->prim;
1.378 + /* determine new bounds of x[j] for down- and up-branches */
1.379 + new_ub = floor(beta);
1.380 + new_lb = ceil(beta);
1.381 + switch (type)
1.382 + { case GLP_FR:
1.383 + dn_type = GLP_UP;
1.384 + up_type = GLP_LO;
1.385 + break;
1.386 + case GLP_LO:
1.387 + xassert(lb <= new_ub);
1.388 + dn_type = (lb == new_ub ? GLP_FX : GLP_DB);
1.389 + xassert(lb + 1.0 <= new_lb);
1.390 + up_type = GLP_LO;
1.391 + break;
1.392 + case GLP_UP:
1.393 + xassert(new_ub <= ub - 1.0);
1.394 + dn_type = GLP_UP;
1.395 + xassert(new_lb <= ub);
1.396 + up_type = (new_lb == ub ? GLP_FX : GLP_DB);
1.397 + break;
1.398 + case GLP_DB:
1.399 + xassert(lb <= new_ub && new_ub <= ub - 1.0);
1.400 + dn_type = (lb == new_ub ? GLP_FX : GLP_DB);
1.401 + xassert(lb + 1.0 <= new_lb && new_lb <= ub);
1.402 + up_type = (new_lb == ub ? GLP_FX : GLP_DB);
1.403 + break;
1.404 + default:
1.405 + xassert(type != type);
1.406 + }
1.407 + /* compute local bounds to LP relaxation for both branches */
1.408 + ios_eval_degrad(T, j, &dn_lp, &up_lp);
1.409 + /* and improve them by rounding */
1.410 + dn_bnd = ios_round_bound(T, dn_lp);
1.411 + up_bnd = ios_round_bound(T, up_lp);
1.412 + /* check local bounds for down- and up-branches */
1.413 + dn_bad = !ios_is_hopeful(T, dn_bnd);
1.414 + up_bad = !ios_is_hopeful(T, up_bnd);
1.415 + if (dn_bad && up_bad)
1.416 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.417 + xprintf("Both down- and up-branches are hopeless\n");
1.418 + ret = 2;
1.419 + goto done;
1.420 + }
1.421 + else if (up_bad)
1.422 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.423 + xprintf("Up-branch is hopeless\n");
1.424 + glp_set_col_bnds(mip, j, dn_type, lb, new_ub);
1.425 + T->curr->lp_obj = dn_lp;
1.426 + if (mip->dir == GLP_MIN)
1.427 + { if (T->curr->bound < dn_bnd)
1.428 + T->curr->bound = dn_bnd;
1.429 + }
1.430 + else if (mip->dir == GLP_MAX)
1.431 + { if (T->curr->bound > dn_bnd)
1.432 + T->curr->bound = dn_bnd;
1.433 + }
1.434 + else
1.435 + xassert(mip != mip);
1.436 + ret = 1;
1.437 + goto done;
1.438 + }
1.439 + else if (dn_bad)
1.440 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.441 + xprintf("Down-branch is hopeless\n");
1.442 + glp_set_col_bnds(mip, j, up_type, new_lb, ub);
1.443 + T->curr->lp_obj = up_lp;
1.444 + if (mip->dir == GLP_MIN)
1.445 + { if (T->curr->bound < up_bnd)
1.446 + T->curr->bound = up_bnd;
1.447 + }
1.448 + else if (mip->dir == GLP_MAX)
1.449 + { if (T->curr->bound > up_bnd)
1.450 + T->curr->bound = up_bnd;
1.451 + }
1.452 + else
1.453 + xassert(mip != mip);
1.454 + ret = 1;
1.455 + goto done;
1.456 + }
1.457 + /* both down- and up-branches seem to be hopeful */
1.458 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.459 + xprintf("Branching on column %d, primal value is %.9e\n",
1.460 + j, beta);
1.461 + /* determine the reference number of the current subproblem */
1.462 + xassert(T->curr != NULL);
1.463 + p = T->curr->p;
1.464 + T->curr->br_var = j;
1.465 + T->curr->br_val = beta;
1.466 + /* freeze the current subproblem */
1.467 + ios_freeze_node(T);
1.468 + /* create two clones of the current subproblem; the first clone
1.469 + begins the down-branch, the second one begins the up-branch */
1.470 + ios_clone_node(T, p, 2, clone);
1.471 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.472 + xprintf("Node %d begins down branch, node %d begins up branch "
1.473 + "\n", clone[1], clone[2]);
1.474 + /* set new upper bound of j-th column in the down-branch */
1.475 + node = T->slot[clone[1]].node;
1.476 + xassert(node != NULL);
1.477 + xassert(node->up != NULL);
1.478 + xassert(node->b_ptr == NULL);
1.479 + node->b_ptr = dmp_get_atom(T->pool, sizeof(IOSBND));
1.480 + node->b_ptr->k = m + j;
1.481 + node->b_ptr->type = (unsigned char)dn_type;
1.482 + node->b_ptr->lb = lb;
1.483 + node->b_ptr->ub = new_ub;
1.484 + node->b_ptr->next = NULL;
1.485 + node->lp_obj = dn_lp;
1.486 + if (mip->dir == GLP_MIN)
1.487 + { if (node->bound < dn_bnd)
1.488 + node->bound = dn_bnd;
1.489 + }
1.490 + else if (mip->dir == GLP_MAX)
1.491 + { if (node->bound > dn_bnd)
1.492 + node->bound = dn_bnd;
1.493 + }
1.494 + else
1.495 + xassert(mip != mip);
1.496 + /* set new lower bound of j-th column in the up-branch */
1.497 + node = T->slot[clone[2]].node;
1.498 + xassert(node != NULL);
1.499 + xassert(node->up != NULL);
1.500 + xassert(node->b_ptr == NULL);
1.501 + node->b_ptr = dmp_get_atom(T->pool, sizeof(IOSBND));
1.502 + node->b_ptr->k = m + j;
1.503 + node->b_ptr->type = (unsigned char)up_type;
1.504 + node->b_ptr->lb = new_lb;
1.505 + node->b_ptr->ub = ub;
1.506 + node->b_ptr->next = NULL;
1.507 + node->lp_obj = up_lp;
1.508 + if (mip->dir == GLP_MIN)
1.509 + { if (node->bound < up_bnd)
1.510 + node->bound = up_bnd;
1.511 + }
1.512 + else if (mip->dir == GLP_MAX)
1.513 + { if (node->bound > up_bnd)
1.514 + node->bound = up_bnd;
1.515 + }
1.516 + else
1.517 + xassert(mip != mip);
1.518 + /* suggest the subproblem to be solved next */
1.519 + xassert(T->child == 0);
1.520 + if (next == GLP_NO_BRNCH)
1.521 + T->child = 0;
1.522 + else if (next == GLP_DN_BRNCH)
1.523 + T->child = clone[1];
1.524 + else if (next == GLP_UP_BRNCH)
1.525 + T->child = clone[2];
1.526 + else
1.527 + xassert(next != next);
1.528 + ret = 0;
1.529 +done: return ret;
1.530 +}
1.531 +
1.532 +/***********************************************************************
1.533 +* cleanup_the_tree - prune hopeless branches from the tree
1.534 +*
1.535 +* This routine walks through the active list and checks the local
1.536 +* bound for every active subproblem. If the local bound indicates that
1.537 +* the subproblem cannot have integer optimal solution better than the
1.538 +* incumbent objective value, the routine deletes such subproblem that,
1.539 +* in turn, involves pruning the corresponding branch of the tree. */
1.540 +
1.541 +static void cleanup_the_tree(glp_tree *T)
1.542 +{ IOSNPD *node, *next_node;
1.543 + int count = 0;
1.544 + /* the global bound must exist */
1.545 + xassert(T->mip->mip_stat == GLP_FEAS);
1.546 + /* walk through the list of active subproblems */
1.547 + for (node = T->head; node != NULL; node = next_node)
1.548 + { /* deleting some active problem node may involve deleting its
1.549 + parents recursively; however, all its parents being created
1.550 + *before* it are always *precede* it in the node list, so
1.551 + the next problem node is never affected by such deletion */
1.552 + next_node = node->next;
1.553 + /* if the branch is hopeless, prune it */
1.554 + if (!is_branch_hopeful(T, node->p))
1.555 + ios_delete_node(T, node->p), count++;
1.556 + }
1.557 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.558 + { if (count == 1)
1.559 + xprintf("One hopeless branch has been pruned\n");
1.560 + else if (count > 1)
1.561 + xprintf("%d hopeless branches have been pruned\n", count);
1.562 + }
1.563 + return;
1.564 +}
1.565 +
1.566 +/**********************************************************************/
1.567 +
1.568 +static void generate_cuts(glp_tree *T)
1.569 +{ /* generate generic cuts with built-in generators */
1.570 + if (!(T->parm->mir_cuts == GLP_ON ||
1.571 + T->parm->gmi_cuts == GLP_ON ||
1.572 + T->parm->cov_cuts == GLP_ON ||
1.573 + T->parm->clq_cuts == GLP_ON)) goto done;
1.574 +#if 1 /* 20/IX-2008 */
1.575 + { int i, max_cuts, added_cuts;
1.576 + max_cuts = T->n;
1.577 + if (max_cuts < 1000) max_cuts = 1000;
1.578 + added_cuts = 0;
1.579 + for (i = T->orig_m+1; i <= T->mip->m; i++)
1.580 + { if (T->mip->row[i]->origin == GLP_RF_CUT)
1.581 + added_cuts++;
1.582 + }
1.583 + /* xprintf("added_cuts = %d\n", added_cuts); */
1.584 + if (added_cuts >= max_cuts) goto done;
1.585 + }
1.586 +#endif
1.587 + /* generate and add to POOL all cuts violated by x* */
1.588 + if (T->parm->gmi_cuts == GLP_ON)
1.589 + { if (T->curr->changed < 5)
1.590 + ios_gmi_gen(T);
1.591 + }
1.592 + if (T->parm->mir_cuts == GLP_ON)
1.593 + { xassert(T->mir_gen != NULL);
1.594 + ios_mir_gen(T, T->mir_gen);
1.595 + }
1.596 + if (T->parm->cov_cuts == GLP_ON)
1.597 + { /* cover cuts works well along with mir cuts */
1.598 + /*if (T->round <= 5)*/
1.599 + ios_cov_gen(T);
1.600 + }
1.601 + if (T->parm->clq_cuts == GLP_ON)
1.602 + { if (T->clq_gen != NULL)
1.603 + { if (T->curr->level == 0 && T->curr->changed < 50 ||
1.604 + T->curr->level > 0 && T->curr->changed < 5)
1.605 + ios_clq_gen(T, T->clq_gen);
1.606 + }
1.607 + }
1.608 +done: return;
1.609 +}
1.610 +
1.611 +/**********************************************************************/
1.612 +
1.613 +static void remove_cuts(glp_tree *T)
1.614 +{ /* remove inactive cuts (some valueable globally valid cut might
1.615 + be saved in the global cut pool) */
1.616 + int i, cnt = 0, *num = NULL;
1.617 + xassert(T->curr != NULL);
1.618 + for (i = T->orig_m+1; i <= T->mip->m; i++)
1.619 + { if (T->mip->row[i]->origin == GLP_RF_CUT &&
1.620 + T->mip->row[i]->level == T->curr->level &&
1.621 + T->mip->row[i]->stat == GLP_BS)
1.622 + { if (num == NULL)
1.623 + num = xcalloc(1+T->mip->m, sizeof(int));
1.624 + num[++cnt] = i;
1.625 + }
1.626 + }
1.627 + if (cnt > 0)
1.628 + { glp_del_rows(T->mip, cnt, num);
1.629 +#if 0
1.630 + xprintf("%d inactive cut(s) removed\n", cnt);
1.631 +#endif
1.632 + xfree(num);
1.633 + xassert(glp_factorize(T->mip) == 0);
1.634 + }
1.635 + return;
1.636 +}
1.637 +
1.638 +/**********************************************************************/
1.639 +
1.640 +static void display_cut_info(glp_tree *T)
1.641 +{ glp_prob *mip = T->mip;
1.642 + int i, gmi = 0, mir = 0, cov = 0, clq = 0, app = 0;
1.643 + for (i = mip->m; i > 0; i--)
1.644 + { GLPROW *row;
1.645 + row = mip->row[i];
1.646 + /* if (row->level < T->curr->level) break; */
1.647 + if (row->origin == GLP_RF_CUT)
1.648 + { if (row->klass == GLP_RF_GMI)
1.649 + gmi++;
1.650 + else if (row->klass == GLP_RF_MIR)
1.651 + mir++;
1.652 + else if (row->klass == GLP_RF_COV)
1.653 + cov++;
1.654 + else if (row->klass == GLP_RF_CLQ)
1.655 + clq++;
1.656 + else
1.657 + app++;
1.658 + }
1.659 + }
1.660 + xassert(T->curr != NULL);
1.661 + if (gmi + mir + cov + clq + app > 0)
1.662 + { xprintf("Cuts on level %d:", T->curr->level);
1.663 + if (gmi > 0) xprintf(" gmi = %d;", gmi);
1.664 + if (mir > 0) xprintf(" mir = %d;", mir);
1.665 + if (cov > 0) xprintf(" cov = %d;", cov);
1.666 + if (clq > 0) xprintf(" clq = %d;", clq);
1.667 + if (app > 0) xprintf(" app = %d;", app);
1.668 + xprintf("\n");
1.669 + }
1.670 + return;
1.671 +}
1.672 +
1.673 +/***********************************************************************
1.674 +* NAME
1.675 +*
1.676 +* ios_driver - branch-and-cut driver
1.677 +*
1.678 +* SYNOPSIS
1.679 +*
1.680 +* #include "glpios.h"
1.681 +* int ios_driver(glp_tree *T);
1.682 +*
1.683 +* DESCRIPTION
1.684 +*
1.685 +* The routine ios_driver is a branch-and-cut driver. It controls the
1.686 +* MIP solution process.
1.687 +*
1.688 +* RETURNS
1.689 +*
1.690 +* 0 The MIP problem instance has been successfully solved. This code
1.691 +* does not necessarily mean that the solver has found optimal
1.692 +* solution. It only means that the solution process was successful.
1.693 +*
1.694 +* GLP_EFAIL
1.695 +* The search was prematurely terminated due to the solver failure.
1.696 +*
1.697 +* GLP_EMIPGAP
1.698 +* The search was prematurely terminated, because the relative mip
1.699 +* gap tolerance has been reached.
1.700 +*
1.701 +* GLP_ETMLIM
1.702 +* The search was prematurely terminated, because the time limit has
1.703 +* been exceeded.
1.704 +*
1.705 +* GLP_ESTOP
1.706 +* The search was prematurely terminated by application. */
1.707 +
1.708 +int ios_driver(glp_tree *T)
1.709 +{ int p, curr_p, p_stat, d_stat, ret;
1.710 +#if 1 /* carry out to glp_tree */
1.711 + int pred_p = 0;
1.712 + /* if the current subproblem has been just created due to
1.713 + branching, pred_p is the reference number of its parent
1.714 + subproblem, otherwise pred_p is zero */
1.715 +#endif
1.716 + glp_long ttt = T->tm_beg;
1.717 +#if 0
1.718 + ((glp_iocp *)T->parm)->msg_lev = GLP_MSG_DBG;
1.719 +#endif
1.720 + /* on entry to the B&B driver it is assumed that the active list
1.721 + contains the only active (i.e. root) subproblem, which is the
1.722 + original MIP problem to be solved */
1.723 +loop: /* main loop starts here */
1.724 + /* at this point the current subproblem does not exist */
1.725 + xassert(T->curr == NULL);
1.726 + /* if the active list is empty, the search is finished */
1.727 + if (T->head == NULL)
1.728 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.729 + xprintf("Active list is empty!\n");
1.730 + xassert(dmp_in_use(T->pool).lo == 0);
1.731 + ret = 0;
1.732 + goto done;
1.733 + }
1.734 + /* select some active subproblem to continue the search */
1.735 + xassert(T->next_p == 0);
1.736 + /* let the application program select subproblem */
1.737 + if (T->parm->cb_func != NULL)
1.738 + { xassert(T->reason == 0);
1.739 + T->reason = GLP_ISELECT;
1.740 + T->parm->cb_func(T, T->parm->cb_info);
1.741 + T->reason = 0;
1.742 + if (T->stop)
1.743 + { ret = GLP_ESTOP;
1.744 + goto done;
1.745 + }
1.746 + }
1.747 + if (T->next_p != 0)
1.748 + { /* the application program has selected something */
1.749 + ;
1.750 + }
1.751 + else if (T->a_cnt == 1)
1.752 + { /* the only active subproblem exists, so select it */
1.753 + xassert(T->head->next == NULL);
1.754 + T->next_p = T->head->p;
1.755 + }
1.756 + else if (T->child != 0)
1.757 + { /* select one of branching childs suggested by the branching
1.758 + heuristic */
1.759 + T->next_p = T->child;
1.760 + }
1.761 + else
1.762 + { /* select active subproblem as specified by the backtracking
1.763 + technique option */
1.764 + T->next_p = ios_choose_node(T);
1.765 + }
1.766 + /* the active subproblem just selected becomes current */
1.767 + ios_revive_node(T, T->next_p);
1.768 + T->next_p = T->child = 0;
1.769 + /* invalidate pred_p, if it is not the reference number of the
1.770 + parent of the current subproblem */
1.771 + if (T->curr->up != NULL && T->curr->up->p != pred_p) pred_p = 0;
1.772 + /* determine the reference number of the current subproblem */
1.773 + p = T->curr->p;
1.774 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.775 + { xprintf("-----------------------------------------------------"
1.776 + "-------------------\n");
1.777 + xprintf("Processing node %d at level %d\n", p, T->curr->level);
1.778 + }
1.779 + /* if it is the root subproblem, initialize cut generators */
1.780 + if (p == 1)
1.781 + { if (T->parm->gmi_cuts == GLP_ON)
1.782 + { if (T->parm->msg_lev >= GLP_MSG_ALL)
1.783 + xprintf("Gomory's cuts enabled\n");
1.784 + }
1.785 + if (T->parm->mir_cuts == GLP_ON)
1.786 + { if (T->parm->msg_lev >= GLP_MSG_ALL)
1.787 + xprintf("MIR cuts enabled\n");
1.788 + xassert(T->mir_gen == NULL);
1.789 + T->mir_gen = ios_mir_init(T);
1.790 + }
1.791 + if (T->parm->cov_cuts == GLP_ON)
1.792 + { if (T->parm->msg_lev >= GLP_MSG_ALL)
1.793 + xprintf("Cover cuts enabled\n");
1.794 + }
1.795 + if (T->parm->clq_cuts == GLP_ON)
1.796 + { xassert(T->clq_gen == NULL);
1.797 + if (T->parm->msg_lev >= GLP_MSG_ALL)
1.798 + xprintf("Clique cuts enabled\n");
1.799 + T->clq_gen = ios_clq_init(T);
1.800 + }
1.801 + }
1.802 +more: /* minor loop starts here */
1.803 + /* at this point the current subproblem needs either to be solved
1.804 + for the first time or re-optimized due to reformulation */
1.805 + /* display current progress of the search */
1.806 + if (T->parm->msg_lev >= GLP_MSG_DBG ||
1.807 + T->parm->msg_lev >= GLP_MSG_ON &&
1.808 + (double)(T->parm->out_frq - 1) <=
1.809 + 1000.0 * xdifftime(xtime(), T->tm_lag))
1.810 + show_progress(T, 0);
1.811 + if (T->parm->msg_lev >= GLP_MSG_ALL &&
1.812 + xdifftime(xtime(), ttt) >= 60.0)
1.813 + { glp_long total;
1.814 + glp_mem_usage(NULL, NULL, &total, NULL);
1.815 + xprintf("Time used: %.1f secs. Memory used: %.1f Mb.\n",
1.816 + xdifftime(xtime(), T->tm_beg), xltod(total) / 1048576.0);
1.817 + ttt = xtime();
1.818 + }
1.819 + /* check the mip gap */
1.820 + if (T->parm->mip_gap > 0.0 &&
1.821 + ios_relative_gap(T) <= T->parm->mip_gap)
1.822 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.823 + xprintf("Relative gap tolerance reached; search terminated "
1.824 + "\n");
1.825 + ret = GLP_EMIPGAP;
1.826 + goto done;
1.827 + }
1.828 + /* check if the time limit has been exhausted */
1.829 + if (T->parm->tm_lim < INT_MAX &&
1.830 + (double)(T->parm->tm_lim - 1) <=
1.831 + 1000.0 * xdifftime(xtime(), T->tm_beg))
1.832 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.833 + xprintf("Time limit exhausted; search terminated\n");
1.834 + ret = GLP_ETMLIM;
1.835 + goto done;
1.836 + }
1.837 + /* let the application program preprocess the subproblem */
1.838 + if (T->parm->cb_func != NULL)
1.839 + { xassert(T->reason == 0);
1.840 + T->reason = GLP_IPREPRO;
1.841 + T->parm->cb_func(T, T->parm->cb_info);
1.842 + T->reason = 0;
1.843 + if (T->stop)
1.844 + { ret = GLP_ESTOP;
1.845 + goto done;
1.846 + }
1.847 + }
1.848 + /* perform basic preprocessing */
1.849 + if (T->parm->pp_tech == GLP_PP_NONE)
1.850 + ;
1.851 + else if (T->parm->pp_tech == GLP_PP_ROOT)
1.852 + { if (T->curr->level == 0)
1.853 + { if (ios_preprocess_node(T, 100))
1.854 + goto fath;
1.855 + }
1.856 + }
1.857 + else if (T->parm->pp_tech == GLP_PP_ALL)
1.858 + { if (ios_preprocess_node(T, T->curr->level == 0 ? 100 : 10))
1.859 + goto fath;
1.860 + }
1.861 + else
1.862 + xassert(T != T);
1.863 + /* preprocessing may improve the global bound */
1.864 + if (!is_branch_hopeful(T, p))
1.865 + { xprintf("*** not tested yet ***\n");
1.866 + goto fath;
1.867 + }
1.868 + /* solve LP relaxation of the current subproblem */
1.869 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.870 + xprintf("Solving LP relaxation...\n");
1.871 + ret = ios_solve_node(T);
1.872 + if (!(ret == 0 || ret == GLP_EOBJLL || ret == GLP_EOBJUL))
1.873 + { if (T->parm->msg_lev >= GLP_MSG_ERR)
1.874 + xprintf("ios_driver: unable to solve current LP relaxation;"
1.875 + " glp_simplex returned %d\n", ret);
1.876 + ret = GLP_EFAIL;
1.877 + goto done;
1.878 + }
1.879 + /* analyze status of the basic solution to LP relaxation found */
1.880 + p_stat = T->mip->pbs_stat;
1.881 + d_stat = T->mip->dbs_stat;
1.882 + if (p_stat == GLP_FEAS && d_stat == GLP_FEAS)
1.883 + { /* LP relaxation has optimal solution */
1.884 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.885 + xprintf("Found optimal solution to LP relaxation\n");
1.886 + }
1.887 + else if (d_stat == GLP_NOFEAS)
1.888 + { /* LP relaxation has no dual feasible solution */
1.889 + /* since the current subproblem cannot have a larger feasible
1.890 + region than its parent, there is something wrong */
1.891 + if (T->parm->msg_lev >= GLP_MSG_ERR)
1.892 + xprintf("ios_driver: current LP relaxation has no dual feas"
1.893 + "ible solution\n");
1.894 + ret = GLP_EFAIL;
1.895 + goto done;
1.896 + }
1.897 + else if (p_stat == GLP_INFEAS && d_stat == GLP_FEAS)
1.898 + { /* LP relaxation has no primal solution which is better than
1.899 + the incumbent objective value */
1.900 + xassert(T->mip->mip_stat == GLP_FEAS);
1.901 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.902 + xprintf("LP relaxation has no solution better than incumben"
1.903 + "t objective value\n");
1.904 + /* prune the branch */
1.905 + goto fath;
1.906 + }
1.907 + else if (p_stat == GLP_NOFEAS)
1.908 + { /* LP relaxation has no primal feasible solution */
1.909 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.910 + xprintf("LP relaxation has no feasible solution\n");
1.911 + /* prune the branch */
1.912 + goto fath;
1.913 + }
1.914 + else
1.915 + { /* other cases cannot appear */
1.916 + xassert(T->mip != T->mip);
1.917 + }
1.918 + /* at this point basic solution to LP relaxation of the current
1.919 + subproblem is optimal */
1.920 + xassert(p_stat == GLP_FEAS && d_stat == GLP_FEAS);
1.921 + xassert(T->curr != NULL);
1.922 + T->curr->lp_obj = T->mip->obj_val;
1.923 + /* thus, it defines a local bound to integer optimal solution of
1.924 + the current subproblem */
1.925 + { double bound = T->mip->obj_val;
1.926 + /* some local bound to the current subproblem could be already
1.927 + set before, so we should only improve it */
1.928 + bound = ios_round_bound(T, bound);
1.929 + if (T->mip->dir == GLP_MIN)
1.930 + { if (T->curr->bound < bound)
1.931 + T->curr->bound = bound;
1.932 + }
1.933 + else if (T->mip->dir == GLP_MAX)
1.934 + { if (T->curr->bound > bound)
1.935 + T->curr->bound = bound;
1.936 + }
1.937 + else
1.938 + xassert(T->mip != T->mip);
1.939 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.940 + xprintf("Local bound is %.9e\n", bound);
1.941 + }
1.942 + /* if the local bound indicates that integer optimal solution of
1.943 + the current subproblem cannot be better than the global bound,
1.944 + prune the branch */
1.945 + if (!is_branch_hopeful(T, p))
1.946 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.947 + xprintf("Current branch is hopeless and can be pruned\n");
1.948 + goto fath;
1.949 + }
1.950 + /* let the application program generate additional rows ("lazy"
1.951 + constraints) */
1.952 + xassert(T->reopt == 0);
1.953 + xassert(T->reinv == 0);
1.954 + if (T->parm->cb_func != NULL)
1.955 + { xassert(T->reason == 0);
1.956 + T->reason = GLP_IROWGEN;
1.957 + T->parm->cb_func(T, T->parm->cb_info);
1.958 + T->reason = 0;
1.959 + if (T->stop)
1.960 + { ret = GLP_ESTOP;
1.961 + goto done;
1.962 + }
1.963 + if (T->reopt)
1.964 + { /* some rows were added; re-optimization is needed */
1.965 + T->reopt = T->reinv = 0;
1.966 + goto more;
1.967 + }
1.968 + if (T->reinv)
1.969 + { /* no rows were added, however, some inactive rows were
1.970 + removed */
1.971 + T->reinv = 0;
1.972 + xassert(glp_factorize(T->mip) == 0);
1.973 + }
1.974 + }
1.975 + /* check if the basic solution is integer feasible */
1.976 + check_integrality(T);
1.977 + /* if the basic solution satisfies to all integrality conditions,
1.978 + it is a new, better integer feasible solution */
1.979 + if (T->curr->ii_cnt == 0)
1.980 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.981 + xprintf("New integer feasible solution found\n");
1.982 + if (T->parm->msg_lev >= GLP_MSG_ALL)
1.983 + display_cut_info(T);
1.984 + record_solution(T);
1.985 + if (T->parm->msg_lev >= GLP_MSG_ON)
1.986 + show_progress(T, 1);
1.987 + /* make the application program happy */
1.988 + if (T->parm->cb_func != NULL)
1.989 + { xassert(T->reason == 0);
1.990 + T->reason = GLP_IBINGO;
1.991 + T->parm->cb_func(T, T->parm->cb_info);
1.992 + T->reason = 0;
1.993 + if (T->stop)
1.994 + { ret = GLP_ESTOP;
1.995 + goto done;
1.996 + }
1.997 + }
1.998 + /* since the current subproblem has been fathomed, prune its
1.999 + branch */
1.1000 + goto fath;
1.1001 + }
1.1002 + /* at this point basic solution to LP relaxation of the current
1.1003 + subproblem is optimal, but integer infeasible */
1.1004 + /* try to fix some non-basic structural variables of integer kind
1.1005 + on their current bounds due to reduced costs */
1.1006 + if (T->mip->mip_stat == GLP_FEAS)
1.1007 + fix_by_red_cost(T);
1.1008 + /* let the application program try to find some solution to the
1.1009 + original MIP with a primal heuristic */
1.1010 + if (T->parm->cb_func != NULL)
1.1011 + { xassert(T->reason == 0);
1.1012 + T->reason = GLP_IHEUR;
1.1013 + T->parm->cb_func(T, T->parm->cb_info);
1.1014 + T->reason = 0;
1.1015 + if (T->stop)
1.1016 + { ret = GLP_ESTOP;
1.1017 + goto done;
1.1018 + }
1.1019 + /* check if the current branch became hopeless */
1.1020 + if (!is_branch_hopeful(T, p))
1.1021 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.1022 + xprintf("Current branch became hopeless and can be prune"
1.1023 + "d\n");
1.1024 + goto fath;
1.1025 + }
1.1026 + }
1.1027 + /* try to find solution with the feasibility pump heuristic */
1.1028 + if (T->parm->fp_heur)
1.1029 + { xassert(T->reason == 0);
1.1030 + T->reason = GLP_IHEUR;
1.1031 + ios_feas_pump(T);
1.1032 + T->reason = 0;
1.1033 + /* check if the current branch became hopeless */
1.1034 + if (!is_branch_hopeful(T, p))
1.1035 + { if (T->parm->msg_lev >= GLP_MSG_DBG)
1.1036 + xprintf("Current branch became hopeless and can be prune"
1.1037 + "d\n");
1.1038 + goto fath;
1.1039 + }
1.1040 + }
1.1041 + /* it's time to generate cutting planes */
1.1042 + xassert(T->local != NULL);
1.1043 + xassert(T->local->size == 0);
1.1044 + /* let the application program generate some cuts; note that it
1.1045 + can add cuts either to the local cut pool or directly to the
1.1046 + current subproblem */
1.1047 + if (T->parm->cb_func != NULL)
1.1048 + { xassert(T->reason == 0);
1.1049 + T->reason = GLP_ICUTGEN;
1.1050 + T->parm->cb_func(T, T->parm->cb_info);
1.1051 + T->reason = 0;
1.1052 + if (T->stop)
1.1053 + { ret = GLP_ESTOP;
1.1054 + goto done;
1.1055 + }
1.1056 + }
1.1057 + /* try to generate generic cuts with built-in generators
1.1058 + (as suggested by Matteo Fischetti et al. the built-in cuts
1.1059 + are not generated at each branching node; an intense attempt
1.1060 + of generating new cuts is only made at the root node, and then
1.1061 + a moderate effort is spent after each backtracking step) */
1.1062 + if (T->curr->level == 0 || pred_p == 0)
1.1063 + { xassert(T->reason == 0);
1.1064 + T->reason = GLP_ICUTGEN;
1.1065 + generate_cuts(T);
1.1066 + T->reason = 0;
1.1067 + }
1.1068 + /* if the local cut pool is not empty, select useful cuts and add
1.1069 + them to the current subproblem */
1.1070 + if (T->local->size > 0)
1.1071 + { xassert(T->reason == 0);
1.1072 + T->reason = GLP_ICUTGEN;
1.1073 + ios_process_cuts(T);
1.1074 + T->reason = 0;
1.1075 + }
1.1076 + /* clear the local cut pool */
1.1077 + ios_clear_pool(T, T->local);
1.1078 + /* perform re-optimization, if necessary */
1.1079 + if (T->reopt)
1.1080 + { T->reopt = 0;
1.1081 + T->curr->changed++;
1.1082 + goto more;
1.1083 + }
1.1084 + /* no cuts were generated; remove inactive cuts */
1.1085 + remove_cuts(T);
1.1086 + if (T->parm->msg_lev >= GLP_MSG_ALL && T->curr->level == 0)
1.1087 + display_cut_info(T);
1.1088 + /* update history information used on pseudocost branching */
1.1089 + if (T->pcost != NULL) ios_pcost_update(T);
1.1090 + /* it's time to perform branching */
1.1091 + xassert(T->br_var == 0);
1.1092 + xassert(T->br_sel == 0);
1.1093 + /* let the application program choose variable to branch on */
1.1094 + if (T->parm->cb_func != NULL)
1.1095 + { xassert(T->reason == 0);
1.1096 + xassert(T->br_var == 0);
1.1097 + xassert(T->br_sel == 0);
1.1098 + T->reason = GLP_IBRANCH;
1.1099 + T->parm->cb_func(T, T->parm->cb_info);
1.1100 + T->reason = 0;
1.1101 + if (T->stop)
1.1102 + { ret = GLP_ESTOP;
1.1103 + goto done;
1.1104 + }
1.1105 + }
1.1106 + /* if nothing has been chosen, choose some variable as specified
1.1107 + by the branching technique option */
1.1108 + if (T->br_var == 0)
1.1109 + T->br_var = ios_choose_var(T, &T->br_sel);
1.1110 + /* perform actual branching */
1.1111 + curr_p = T->curr->p;
1.1112 + ret = branch_on(T, T->br_var, T->br_sel);
1.1113 + T->br_var = T->br_sel = 0;
1.1114 + if (ret == 0)
1.1115 + { /* both branches have been created */
1.1116 + pred_p = curr_p;
1.1117 + goto loop;
1.1118 + }
1.1119 + else if (ret == 1)
1.1120 + { /* one branch is hopeless and has been pruned, so now the
1.1121 + current subproblem is other branch */
1.1122 + /* the current subproblem should be considered as a new one,
1.1123 + since one bound of the branching variable was changed */
1.1124 + T->curr->solved = T->curr->changed = 0;
1.1125 + goto more;
1.1126 + }
1.1127 + else if (ret == 2)
1.1128 + { /* both branches are hopeless and have been pruned; new
1.1129 + subproblem selection is needed to continue the search */
1.1130 + goto fath;
1.1131 + }
1.1132 + else
1.1133 + xassert(ret != ret);
1.1134 +fath: /* the current subproblem has been fathomed */
1.1135 + if (T->parm->msg_lev >= GLP_MSG_DBG)
1.1136 + xprintf("Node %d fathomed\n", p);
1.1137 + /* freeze the current subproblem */
1.1138 + ios_freeze_node(T);
1.1139 + /* and prune the corresponding branch of the tree */
1.1140 + ios_delete_node(T, p);
1.1141 + /* if a new integer feasible solution has just been found, other
1.1142 + branches may become hopeless and therefore must be pruned */
1.1143 + if (T->mip->mip_stat == GLP_FEAS) cleanup_the_tree(T);
1.1144 + /* new subproblem selection is needed due to backtracking */
1.1145 + pred_p = 0;
1.1146 + goto loop;
1.1147 +done: /* display progress of the search on exit from the solver */
1.1148 + if (T->parm->msg_lev >= GLP_MSG_ON)
1.1149 + show_progress(T, 0);
1.1150 + if (T->mir_gen != NULL)
1.1151 + ios_mir_term(T->mir_gen), T->mir_gen = NULL;
1.1152 + if (T->clq_gen != NULL)
1.1153 + ios_clq_term(T->clq_gen), T->clq_gen = NULL;
1.1154 + /* return to the calling program */
1.1155 + return ret;
1.1156 +}
1.1157 +
1.1158 +/* eof */