1 /* glpapi09.c (mixed integer programming routines) */
3 /***********************************************************************
4 * This code is part of GLPK (GNU Linear Programming Kit).
6 * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
7 * 2009, 2010 Andrew Makhorin, Department for Applied Informatics,
8 * Moscow Aviation Institute, Moscow, Russia. All rights reserved.
9 * E-mail: <mao@gnu.org>.
11 * GLPK is free software: you can redistribute it and/or modify it
12 * under the terms of the GNU General Public License as published by
13 * the Free Software Foundation, either version 3 of the License, or
14 * (at your option) any later version.
16 * GLPK is distributed in the hope that it will be useful, but WITHOUT
17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
19 * License for more details.
21 * You should have received a copy of the GNU General Public License
22 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.
23 ***********************************************************************/
28 /***********************************************************************
31 * glp_set_col_kind - set (change) column kind
35 * void glp_set_col_kind(glp_prob *mip, int j, int kind);
39 * The routine glp_set_col_kind sets (changes) the kind of j-th column
40 * (structural variable) as specified by the parameter kind:
42 * GLP_CV - continuous variable;
43 * GLP_IV - integer variable;
44 * GLP_BV - binary variable. */
46 void glp_set_col_kind(glp_prob *mip, int j, int kind)
48 if (!(1 <= j && j <= mip->n))
49 xerror("glp_set_col_kind: j = %d; column number out of range\n"
61 if (!(col->type == GLP_DB && col->lb == 0.0 && col->ub ==
62 1.0)) glp_set_col_bnds(mip, j, GLP_DB, 0.0, 1.0);
65 xerror("glp_set_col_kind: j = %d; kind = %d; invalid column"
71 /***********************************************************************
74 * glp_get_col_kind - retrieve column kind
78 * int glp_get_col_kind(glp_prob *mip, int j);
82 * The routine glp_get_col_kind returns the kind of j-th column, i.e.
83 * the kind of corresponding structural variable, as follows:
85 * GLP_CV - continuous variable;
86 * GLP_IV - integer variable;
87 * GLP_BV - binary variable */
89 int glp_get_col_kind(glp_prob *mip, int j)
92 if (!(1 <= j && j <= mip->n))
93 xerror("glp_get_col_kind: j = %d; column number out of range\n"
101 if (col->type == GLP_DB && col->lb == 0.0 && col->ub == 1.0)
105 xassert(kind != kind);
110 /***********************************************************************
113 * glp_get_num_int - retrieve number of integer columns
117 * int glp_get_num_int(glp_prob *mip);
121 * The routine glp_get_num_int returns the current number of columns,
122 * which are marked as integer. */
124 int glp_get_num_int(glp_prob *mip)
127 for (j = 1; j <= mip->n; j++)
129 if (col->kind == GLP_IV) count++;
134 /***********************************************************************
137 * glp_get_num_bin - retrieve number of binary columns
141 * int glp_get_num_bin(glp_prob *mip);
145 * The routine glp_get_num_bin returns the current number of columns,
146 * which are marked as binary. */
148 int glp_get_num_bin(glp_prob *mip)
151 for (j = 1; j <= mip->n; j++)
153 if (col->kind == GLP_IV && col->type == GLP_DB && col->lb ==
154 0.0 && col->ub == 1.0) count++;
159 /***********************************************************************
162 * glp_intopt - solve MIP problem with the branch-and-bound method
166 * int glp_intopt(glp_prob *P, const glp_iocp *parm);
170 * The routine glp_intopt is a driver to the MIP solver based on the
171 * branch-and-bound method.
173 * On entry the problem object should contain optimal solution to LP
174 * relaxation (which can be obtained with the routine glp_simplex).
176 * The MIP solver has a set of control parameters. Values of the control
177 * parameters can be passed in a structure glp_iocp, which the parameter
180 * The parameter parm can be specified as NULL, in which case the MIP
181 * solver uses default settings.
185 * 0 The MIP problem instance has been successfully solved. This code
186 * does not necessarily mean that the solver has found optimal
187 * solution. It only means that the solution process was successful.
190 * Unable to start the search, because some double-bounded variables
191 * have incorrect bounds or some integer variables have non-integer
192 * (fractional) bounds.
195 * Unable to start the search, because optimal basis for initial LP
196 * relaxation is not provided.
199 * The search was prematurely terminated due to the solver failure.
202 * The search was prematurely terminated, because the relative mip
203 * gap tolerance has been reached.
206 * The search was prematurely terminated, because the time limit has
210 * The MIP problem instance has no primal feasible solution (only if
211 * the MIP presolver is used).
214 * LP relaxation of the MIP problem instance has no dual feasible
215 * solution (only if the MIP presolver is used).
218 * The search was prematurely terminated by application. */
220 static int solve_mip(glp_prob *P, const glp_iocp *parm)
221 { /* solve MIP directly without using the preprocessor */
224 /* optimal basis to LP relaxation must be provided */
225 if (glp_get_status(P) != GLP_OPT)
226 { if (parm->msg_lev >= GLP_MSG_ERR)
227 xprintf("glp_intopt: optimal basis to initial LP relaxation"
232 /* it seems all is ok */
233 if (parm->msg_lev >= GLP_MSG_ALL)
234 xprintf("Integer optimization begins...\n");
235 /* create the branch-and-bound tree */
236 T = ios_create_tree(P, parm);
237 /* solve the problem instance */
239 /* delete the branch-and-bound tree */
241 /* analyze exit code reported by the mip driver */
243 { if (P->mip_stat == GLP_FEAS)
244 { if (parm->msg_lev >= GLP_MSG_ALL)
245 xprintf("INTEGER OPTIMAL SOLUTION FOUND\n");
246 P->mip_stat = GLP_OPT;
249 { if (parm->msg_lev >= GLP_MSG_ALL)
250 xprintf("PROBLEM HAS NO INTEGER FEASIBLE SOLUTION\n");
251 P->mip_stat = GLP_NOFEAS;
254 else if (ret == GLP_EMIPGAP)
255 { if (parm->msg_lev >= GLP_MSG_ALL)
256 xprintf("RELATIVE MIP GAP TOLERANCE REACHED; SEARCH TERMINA"
259 else if (ret == GLP_ETMLIM)
260 { if (parm->msg_lev >= GLP_MSG_ALL)
261 xprintf("TIME LIMIT EXCEEDED; SEARCH TERMINATED\n");
263 else if (ret == GLP_EFAIL)
264 { if (parm->msg_lev >= GLP_MSG_ERR)
265 xprintf("glp_intopt: cannot solve current LP relaxation\n");
267 else if (ret == GLP_ESTOP)
268 { if (parm->msg_lev >= GLP_MSG_ALL)
269 xprintf("SEARCH TERMINATED BY APPLICATION\n");
276 static int preprocess_and_solve_mip(glp_prob *P, const glp_iocp *parm)
277 { /* solve MIP using the preprocessor */
278 ENV *env = get_env_ptr();
279 int term_out = env->term_out;
281 glp_prob *mip = NULL;
285 if (parm->msg_lev >= GLP_MSG_ALL)
286 xprintf("Preprocessing...\n");
287 /* create preprocessor workspace */
288 npp = npp_create_wksp();
289 /* load original problem into the preprocessor workspace */
290 npp_load_prob(npp, P, GLP_OFF, GLP_MIP, GLP_OFF);
291 /* process MIP prior to applying the branch-and-bound method */
292 if (!term_out || parm->msg_lev < GLP_MSG_ALL)
293 env->term_out = GLP_OFF;
295 env->term_out = GLP_ON;
296 ret = npp_integer(npp, parm);
297 env->term_out = term_out;
300 else if (ret == GLP_ENOPFS)
301 { if (parm->msg_lev >= GLP_MSG_ALL)
302 xprintf("PROBLEM HAS NO PRIMAL FEASIBLE SOLUTION\n");
304 else if (ret == GLP_ENODFS)
305 { if (parm->msg_lev >= GLP_MSG_ALL)
306 xprintf("LP RELAXATION HAS NO DUAL FEASIBLE SOLUTION\n");
310 if (ret != 0) goto done;
311 /* build transformed MIP */
312 mip = glp_create_prob();
313 npp_build_prob(npp, mip);
314 /* if the transformed MIP is empty, it has empty solution, which
316 if (mip->m == 0 && mip->n == 0)
317 { mip->mip_stat = GLP_OPT;
318 mip->mip_obj = mip->c0;
319 if (parm->msg_lev >= GLP_MSG_ALL)
320 { xprintf("Objective value = %17.9e\n", mip->mip_obj);
321 xprintf("INTEGER OPTIMAL SOLUTION FOUND BY MIP PREPROCESSOR"
326 /* display some statistics */
327 if (parm->msg_lev >= GLP_MSG_ALL)
328 { int ni = glp_get_num_int(mip);
329 int nb = glp_get_num_bin(mip);
331 xprintf("%d row%s, %d column%s, %d non-zero%s\n",
332 mip->m, mip->m == 1 ? "" : "s", mip->n, mip->n == 1 ? "" :
333 "s", mip->nnz, mip->nnz == 1 ? "" : "s");
335 strcpy(s, "none of");
336 else if (ni == 1 && nb == 1)
343 sprintf(s, "%d of", nb);
344 xprintf("%d integer variable%s, %s which %s binary\n",
345 ni, ni == 1 ? "" : "s", s, nb == 1 ? "is" : "are");
347 /* inherit basis factorization control parameters */
348 glp_get_bfcp(P, &bfcp);
349 glp_set_bfcp(mip, &bfcp);
350 /* scale the transformed problem */
351 if (!term_out || parm->msg_lev < GLP_MSG_ALL)
352 env->term_out = GLP_OFF;
354 env->term_out = GLP_ON;
356 GLP_SF_GM | GLP_SF_EQ | GLP_SF_2N | GLP_SF_SKIP);
357 env->term_out = term_out;
358 /* build advanced initial basis */
359 if (!term_out || parm->msg_lev < GLP_MSG_ALL)
360 env->term_out = GLP_OFF;
362 env->term_out = GLP_ON;
363 glp_adv_basis(mip, 0);
364 env->term_out = term_out;
365 /* solve initial LP relaxation */
366 if (parm->msg_lev >= GLP_MSG_ALL)
367 xprintf("Solving LP relaxation...\n");
368 glp_init_smcp(&smcp);
369 smcp.msg_lev = parm->msg_lev;
370 mip->it_cnt = P->it_cnt;
371 ret = glp_simplex(mip, &smcp);
372 P->it_cnt = mip->it_cnt;
374 { if (parm->msg_lev >= GLP_MSG_ERR)
375 xprintf("glp_intopt: cannot solve LP relaxation\n");
379 /* check status of the basic solution */
380 ret = glp_get_status(mip);
383 else if (ret == GLP_NOFEAS)
385 else if (ret == GLP_UNBND)
389 if (ret != 0) goto done;
390 /* solve the transformed MIP */
391 mip->it_cnt = P->it_cnt;
392 ret = solve_mip(mip, parm);
393 P->it_cnt = mip->it_cnt;
394 /* only integer feasible solution can be postprocessed */
395 if (!(mip->mip_stat == GLP_OPT || mip->mip_stat == GLP_FEAS))
396 { P->mip_stat = mip->mip_stat;
399 /* postprocess solution from the transformed MIP */
400 post: npp_postprocess(npp, mip);
401 /* the transformed MIP is no longer needed */
402 glp_delete_prob(mip), mip = NULL;
403 /* store solution to the original problem */
404 npp_unload_sol(npp, P);
405 done: /* delete the transformed MIP, if it exists */
406 if (mip != NULL) glp_delete_prob(mip);
407 /* delete preprocessor workspace */
408 npp_delete_wksp(npp);
412 #ifndef HAVE_ALIEN_SOLVER /* 28/V-2010 */
413 int _glp_intopt1(glp_prob *P, const glp_iocp *parm)
415 xassert(parm == parm);
416 xprintf("glp_intopt: no alien solver is available\n");
421 int glp_intopt(glp_prob *P, const glp_iocp *parm)
422 { /* solve MIP problem with the branch-and-bound method */
425 /* check problem object */
426 if (P == NULL || P->magic != GLP_PROB_MAGIC)
427 xerror("glp_intopt: P = %p; invalid problem object\n", P);
429 xerror("glp_intopt: operation not allowed\n");
430 /* check control parameters */
432 parm = &_parm, glp_init_iocp((glp_iocp *)parm);
433 if (!(parm->msg_lev == GLP_MSG_OFF ||
434 parm->msg_lev == GLP_MSG_ERR ||
435 parm->msg_lev == GLP_MSG_ON ||
436 parm->msg_lev == GLP_MSG_ALL ||
437 parm->msg_lev == GLP_MSG_DBG))
438 xerror("glp_intopt: msg_lev = %d; invalid parameter\n",
440 if (!(parm->br_tech == GLP_BR_FFV ||
441 parm->br_tech == GLP_BR_LFV ||
442 parm->br_tech == GLP_BR_MFV ||
443 parm->br_tech == GLP_BR_DTH ||
444 parm->br_tech == GLP_BR_PCH))
445 xerror("glp_intopt: br_tech = %d; invalid parameter\n",
447 if (!(parm->bt_tech == GLP_BT_DFS ||
448 parm->bt_tech == GLP_BT_BFS ||
449 parm->bt_tech == GLP_BT_BLB ||
450 parm->bt_tech == GLP_BT_BPH))
451 xerror("glp_intopt: bt_tech = %d; invalid parameter\n",
453 if (!(0.0 < parm->tol_int && parm->tol_int < 1.0))
454 xerror("glp_intopt: tol_int = %g; invalid parameter\n",
456 if (!(0.0 < parm->tol_obj && parm->tol_obj < 1.0))
457 xerror("glp_intopt: tol_obj = %g; invalid parameter\n",
459 if (parm->tm_lim < 0)
460 xerror("glp_intopt: tm_lim = %d; invalid parameter\n",
462 if (parm->out_frq < 0)
463 xerror("glp_intopt: out_frq = %d; invalid parameter\n",
465 if (parm->out_dly < 0)
466 xerror("glp_intopt: out_dly = %d; invalid parameter\n",
468 if (!(0 <= parm->cb_size && parm->cb_size <= 256))
469 xerror("glp_intopt: cb_size = %d; invalid parameter\n",
471 if (!(parm->pp_tech == GLP_PP_NONE ||
472 parm->pp_tech == GLP_PP_ROOT ||
473 parm->pp_tech == GLP_PP_ALL))
474 xerror("glp_intopt: pp_tech = %d; invalid parameter\n",
476 if (parm->mip_gap < 0.0)
477 xerror("glp_intopt: mip_gap = %g; invalid parameter\n",
479 if (!(parm->mir_cuts == GLP_ON || parm->mir_cuts == GLP_OFF))
480 xerror("glp_intopt: mir_cuts = %d; invalid parameter\n",
482 if (!(parm->gmi_cuts == GLP_ON || parm->gmi_cuts == GLP_OFF))
483 xerror("glp_intopt: gmi_cuts = %d; invalid parameter\n",
485 if (!(parm->cov_cuts == GLP_ON || parm->cov_cuts == GLP_OFF))
486 xerror("glp_intopt: cov_cuts = %d; invalid parameter\n",
488 if (!(parm->clq_cuts == GLP_ON || parm->clq_cuts == GLP_OFF))
489 xerror("glp_intopt: clq_cuts = %d; invalid parameter\n",
491 if (!(parm->presolve == GLP_ON || parm->presolve == GLP_OFF))
492 xerror("glp_intopt: presolve = %d; invalid parameter\n",
494 if (!(parm->binarize == GLP_ON || parm->binarize == GLP_OFF))
495 xerror("glp_intopt: binarize = %d; invalid parameter\n",
497 if (!(parm->fp_heur == GLP_ON || parm->fp_heur == GLP_OFF))
498 xerror("glp_intopt: fp_heur = %d; invalid parameter\n",
500 #if 1 /* 28/V-2010 */
501 if (!(parm->alien == GLP_ON || parm->alien == GLP_OFF))
502 xerror("glp_intopt: alien = %d; invalid parameter\n",
505 /* integer solution is currently undefined */
506 P->mip_stat = GLP_UNDEF;
508 /* check bounds of double-bounded variables */
509 for (i = 1; i <= P->m; i++)
510 { GLPROW *row = P->row[i];
511 if (row->type == GLP_DB && row->lb >= row->ub)
512 { if (parm->msg_lev >= GLP_MSG_ERR)
513 xprintf("glp_intopt: row %d: lb = %g, ub = %g; incorrect"
514 " bounds\n", i, row->lb, row->ub);
519 for (j = 1; j <= P->n; j++)
520 { GLPCOL *col = P->col[j];
521 if (col->type == GLP_DB && col->lb >= col->ub)
522 { if (parm->msg_lev >= GLP_MSG_ERR)
523 xprintf("glp_intopt: column %d: lb = %g, ub = %g; incorr"
524 "ect bounds\n", j, col->lb, col->ub);
529 /* bounds of all integer variables must be integral */
530 for (j = 1; j <= P->n; j++)
531 { GLPCOL *col = P->col[j];
532 if (col->kind != GLP_IV) continue;
533 if (col->type == GLP_LO || col->type == GLP_DB)
534 { if (col->lb != floor(col->lb))
535 { if (parm->msg_lev >= GLP_MSG_ERR)
536 xprintf("glp_intopt: integer column %d has non-intege"
537 "r lower bound %g\n", j, col->lb);
542 if (col->type == GLP_UP || col->type == GLP_DB)
543 { if (col->ub != floor(col->ub))
544 { if (parm->msg_lev >= GLP_MSG_ERR)
545 xprintf("glp_intopt: integer column %d has non-intege"
546 "r upper bound %g\n", j, col->ub);
551 if (col->type == GLP_FX)
552 { if (col->lb != floor(col->lb))
553 { if (parm->msg_lev >= GLP_MSG_ERR)
554 xprintf("glp_intopt: integer column %d has non-intege"
555 "r fixed value %g\n", j, col->lb);
561 /* solve MIP problem */
562 if (parm->msg_lev >= GLP_MSG_ALL)
563 { int ni = glp_get_num_int(P);
564 int nb = glp_get_num_bin(P);
566 xprintf("GLPK Integer Optimizer, v%s\n", glp_version());
567 xprintf("%d row%s, %d column%s, %d non-zero%s\n",
568 P->m, P->m == 1 ? "" : "s", P->n, P->n == 1 ? "" : "s",
569 P->nnz, P->nnz == 1 ? "" : "s");
571 strcpy(s, "none of");
572 else if (ni == 1 && nb == 1)
579 sprintf(s, "%d of", nb);
580 xprintf("%d integer variable%s, %s which %s binary\n",
581 ni, ni == 1 ? "" : "s", s, nb == 1 ? "is" : "are");
583 #if 1 /* 28/V-2010 */
585 { /* use alien integer optimizer */
586 ret = _glp_intopt1(P, parm);
591 ret = solve_mip(P, parm);
593 ret = preprocess_and_solve_mip(P, parm);
594 done: /* return to the application program */
598 /***********************************************************************
601 * glp_init_iocp - initialize integer optimizer control parameters
605 * void glp_init_iocp(glp_iocp *parm);
609 * The routine glp_init_iocp initializes control parameters, which are
610 * used by the integer optimizer, with default values.
612 * Default values of the control parameters are stored in a glp_iocp
613 * structure, which the parameter parm points to. */
615 void glp_init_iocp(glp_iocp *parm)
616 { parm->msg_lev = GLP_MSG_ALL;
617 parm->br_tech = GLP_BR_DTH;
618 parm->bt_tech = GLP_BT_BLB;
619 parm->tol_int = 1e-5;
620 parm->tol_obj = 1e-7;
621 parm->tm_lim = INT_MAX;
622 parm->out_frq = 5000;
623 parm->out_dly = 10000;
624 parm->cb_func = NULL;
625 parm->cb_info = NULL;
627 parm->pp_tech = GLP_PP_ALL;
629 parm->mir_cuts = GLP_OFF;
630 parm->gmi_cuts = GLP_OFF;
631 parm->cov_cuts = GLP_OFF;
632 parm->clq_cuts = GLP_OFF;
633 parm->presolve = GLP_OFF;
634 parm->binarize = GLP_OFF;
635 parm->fp_heur = GLP_OFF;
636 #if 1 /* 28/V-2010 */
637 parm->alien = GLP_OFF;
642 /***********************************************************************
645 * glp_mip_status - retrieve status of MIP solution
649 * int glp_mip_status(glp_prob *mip);
653 * The routine lpx_mip_status reports the status of MIP solution found
654 * by the branch-and-bound solver as follows:
656 * GLP_UNDEF - MIP solution is undefined;
657 * GLP_OPT - MIP solution is integer optimal;
658 * GLP_FEAS - MIP solution is integer feasible but its optimality
659 * (or non-optimality) has not been proven, perhaps due to
660 * premature termination of the search;
661 * GLP_NOFEAS - problem has no integer feasible solution (proven by the
664 int glp_mip_status(glp_prob *mip)
665 { int mip_stat = mip->mip_stat;
669 /***********************************************************************
672 * glp_mip_obj_val - retrieve objective value (MIP solution)
676 * double glp_mip_obj_val(glp_prob *mip);
680 * The routine glp_mip_obj_val returns value of the objective function
681 * for MIP solution. */
683 double glp_mip_obj_val(glp_prob *mip)
684 { /*struct LPXCPS *cps = mip->cps;*/
687 /*if (cps->round && fabs(z) < 1e-9) z = 0.0;*/
691 /***********************************************************************
694 * glp_mip_row_val - retrieve row value (MIP solution)
698 * double glp_mip_row_val(glp_prob *mip, int i);
702 * The routine glp_mip_row_val returns value of the auxiliary variable
703 * associated with i-th row. */
705 double glp_mip_row_val(glp_prob *mip, int i)
706 { /*struct LPXCPS *cps = mip->cps;*/
708 if (!(1 <= i && i <= mip->m))
709 xerror("glp_mip_row_val: i = %d; row number out of range\n", i)
711 mipx = mip->row[i]->mipx;
712 /*if (cps->round && fabs(mipx) < 1e-9) mipx = 0.0;*/
716 /***********************************************************************
719 * glp_mip_col_val - retrieve column value (MIP solution)
723 * double glp_mip_col_val(glp_prob *mip, int j);
727 * The routine glp_mip_col_val returns value of the structural variable
728 * associated with j-th column. */
730 double glp_mip_col_val(glp_prob *mip, int j)
731 { /*struct LPXCPS *cps = mip->cps;*/
733 if (!(1 <= j && j <= mip->n))
734 xerror("glp_mip_col_val: j = %d; column number out of range\n",
736 mipx = mip->col[j]->mipx;
737 /*if (cps->round && fabs(mipx) < 1e-9) mipx = 0.0;*/