1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/glpapi06.c Mon Dec 06 13:09:21 2010 +0100
1.3 @@ -0,0 +1,806 @@
1.4 +/* glpapi06.c (simplex method routines) */
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 +#include "glpnpp.h"
1.30 +#include "glpspx.h"
1.31 +
1.32 +/***********************************************************************
1.33 +* NAME
1.34 +*
1.35 +* glp_simplex - solve LP problem with the simplex method
1.36 +*
1.37 +* SYNOPSIS
1.38 +*
1.39 +* int glp_simplex(glp_prob *P, const glp_smcp *parm);
1.40 +*
1.41 +* DESCRIPTION
1.42 +*
1.43 +* The routine glp_simplex is a driver to the LP solver based on the
1.44 +* simplex method. This routine retrieves problem data from the
1.45 +* specified problem object, calls the solver to solve the problem
1.46 +* instance, and stores results of computations back into the problem
1.47 +* object.
1.48 +*
1.49 +* The simplex solver has a set of control parameters. Values of the
1.50 +* control parameters can be passed in a structure glp_smcp, which the
1.51 +* parameter parm points to.
1.52 +*
1.53 +* The parameter parm can be specified as NULL, in which case the LP
1.54 +* solver uses default settings.
1.55 +*
1.56 +* RETURNS
1.57 +*
1.58 +* 0 The LP problem instance has been successfully solved. This code
1.59 +* does not necessarily mean that the solver has found optimal
1.60 +* solution. It only means that the solution process was successful.
1.61 +*
1.62 +* GLP_EBADB
1.63 +* Unable to start the search, because the initial basis specified
1.64 +* in the problem object is invalid--the number of basic (auxiliary
1.65 +* and structural) variables is not the same as the number of rows in
1.66 +* the problem object.
1.67 +*
1.68 +* GLP_ESING
1.69 +* Unable to start the search, because the basis matrix correspodning
1.70 +* to the initial basis is singular within the working precision.
1.71 +*
1.72 +* GLP_ECOND
1.73 +* Unable to start the search, because the basis matrix correspodning
1.74 +* to the initial basis is ill-conditioned, i.e. its condition number
1.75 +* is too large.
1.76 +*
1.77 +* GLP_EBOUND
1.78 +* Unable to start the search, because some double-bounded variables
1.79 +* have incorrect bounds.
1.80 +*
1.81 +* GLP_EFAIL
1.82 +* The search was prematurely terminated due to the solver failure.
1.83 +*
1.84 +* GLP_EOBJLL
1.85 +* The search was prematurely terminated, because the objective
1.86 +* function being maximized has reached its lower limit and continues
1.87 +* decreasing (dual simplex only).
1.88 +*
1.89 +* GLP_EOBJUL
1.90 +* The search was prematurely terminated, because the objective
1.91 +* function being minimized has reached its upper limit and continues
1.92 +* increasing (dual simplex only).
1.93 +*
1.94 +* GLP_EITLIM
1.95 +* The search was prematurely terminated, because the simplex
1.96 +* iteration limit has been exceeded.
1.97 +*
1.98 +* GLP_ETMLIM
1.99 +* The search was prematurely terminated, because the time limit has
1.100 +* been exceeded.
1.101 +*
1.102 +* GLP_ENOPFS
1.103 +* The LP problem instance has no primal feasible solution (only if
1.104 +* the LP presolver is used).
1.105 +*
1.106 +* GLP_ENODFS
1.107 +* The LP problem instance has no dual feasible solution (only if the
1.108 +* LP presolver is used). */
1.109 +
1.110 +static void trivial_lp(glp_prob *P, const glp_smcp *parm)
1.111 +{ /* solve trivial LP which has empty constraint matrix */
1.112 + GLPROW *row;
1.113 + GLPCOL *col;
1.114 + int i, j;
1.115 + double p_infeas, d_infeas, zeta;
1.116 + P->valid = 0;
1.117 + P->pbs_stat = P->dbs_stat = GLP_FEAS;
1.118 + P->obj_val = P->c0;
1.119 + P->some = 0;
1.120 + p_infeas = d_infeas = 0.0;
1.121 + /* make all auxiliary variables basic */
1.122 + for (i = 1; i <= P->m; i++)
1.123 + { row = P->row[i];
1.124 + row->stat = GLP_BS;
1.125 + row->prim = row->dual = 0.0;
1.126 + /* check primal feasibility */
1.127 + if (row->type == GLP_LO || row->type == GLP_DB ||
1.128 + row->type == GLP_FX)
1.129 + { /* row has lower bound */
1.130 + if (row->lb > + parm->tol_bnd)
1.131 + { P->pbs_stat = GLP_NOFEAS;
1.132 + if (P->some == 0 && parm->meth != GLP_PRIMAL)
1.133 + P->some = i;
1.134 + }
1.135 + if (p_infeas < + row->lb)
1.136 + p_infeas = + row->lb;
1.137 + }
1.138 + if (row->type == GLP_UP || row->type == GLP_DB ||
1.139 + row->type == GLP_FX)
1.140 + { /* row has upper bound */
1.141 + if (row->ub < - parm->tol_bnd)
1.142 + { P->pbs_stat = GLP_NOFEAS;
1.143 + if (P->some == 0 && parm->meth != GLP_PRIMAL)
1.144 + P->some = i;
1.145 + }
1.146 + if (p_infeas < - row->ub)
1.147 + p_infeas = - row->ub;
1.148 + }
1.149 + }
1.150 + /* determine scale factor for the objective row */
1.151 + zeta = 1.0;
1.152 + for (j = 1; j <= P->n; j++)
1.153 + { col = P->col[j];
1.154 + if (zeta < fabs(col->coef)) zeta = fabs(col->coef);
1.155 + }
1.156 + zeta = (P->dir == GLP_MIN ? +1.0 : -1.0) / zeta;
1.157 + /* make all structural variables non-basic */
1.158 + for (j = 1; j <= P->n; j++)
1.159 + { col = P->col[j];
1.160 + if (col->type == GLP_FR)
1.161 + col->stat = GLP_NF, col->prim = 0.0;
1.162 + else if (col->type == GLP_LO)
1.163 +lo: col->stat = GLP_NL, col->prim = col->lb;
1.164 + else if (col->type == GLP_UP)
1.165 +up: col->stat = GLP_NU, col->prim = col->ub;
1.166 + else if (col->type == GLP_DB)
1.167 + { if (zeta * col->coef > 0.0)
1.168 + goto lo;
1.169 + else if (zeta * col->coef < 0.0)
1.170 + goto up;
1.171 + else if (fabs(col->lb) <= fabs(col->ub))
1.172 + goto lo;
1.173 + else
1.174 + goto up;
1.175 + }
1.176 + else if (col->type == GLP_FX)
1.177 + col->stat = GLP_NS, col->prim = col->lb;
1.178 + col->dual = col->coef;
1.179 + P->obj_val += col->coef * col->prim;
1.180 + /* check dual feasibility */
1.181 + if (col->type == GLP_FR || col->type == GLP_LO)
1.182 + { /* column has no upper bound */
1.183 + if (zeta * col->dual < - parm->tol_dj)
1.184 + { P->dbs_stat = GLP_NOFEAS;
1.185 + if (P->some == 0 && parm->meth == GLP_PRIMAL)
1.186 + P->some = P->m + j;
1.187 + }
1.188 + if (d_infeas < - zeta * col->dual)
1.189 + d_infeas = - zeta * col->dual;
1.190 + }
1.191 + if (col->type == GLP_FR || col->type == GLP_UP)
1.192 + { /* column has no lower bound */
1.193 + if (zeta * col->dual > + parm->tol_dj)
1.194 + { P->dbs_stat = GLP_NOFEAS;
1.195 + if (P->some == 0 && parm->meth == GLP_PRIMAL)
1.196 + P->some = P->m + j;
1.197 + }
1.198 + if (d_infeas < + zeta * col->dual)
1.199 + d_infeas = + zeta * col->dual;
1.200 + }
1.201 + }
1.202 + /* simulate the simplex solver output */
1.203 + if (parm->msg_lev >= GLP_MSG_ON && parm->out_dly == 0)
1.204 + { xprintf("~%6d: obj = %17.9e infeas = %10.3e\n", P->it_cnt,
1.205 + P->obj_val, parm->meth == GLP_PRIMAL ? p_infeas : d_infeas);
1.206 + }
1.207 + if (parm->msg_lev >= GLP_MSG_ALL && parm->out_dly == 0)
1.208 + { if (P->pbs_stat == GLP_FEAS && P->dbs_stat == GLP_FEAS)
1.209 + xprintf("OPTIMAL SOLUTION FOUND\n");
1.210 + else if (P->pbs_stat == GLP_NOFEAS)
1.211 + xprintf("PROBLEM HAS NO FEASIBLE SOLUTION\n");
1.212 + else if (parm->meth == GLP_PRIMAL)
1.213 + xprintf("PROBLEM HAS UNBOUNDED SOLUTION\n");
1.214 + else
1.215 + xprintf("PROBLEM HAS NO DUAL FEASIBLE SOLUTION\n");
1.216 + }
1.217 + return;
1.218 +}
1.219 +
1.220 +static int solve_lp(glp_prob *P, const glp_smcp *parm)
1.221 +{ /* solve LP directly without using the preprocessor */
1.222 + int ret;
1.223 + if (!glp_bf_exists(P))
1.224 + { ret = glp_factorize(P);
1.225 + if (ret == 0)
1.226 + ;
1.227 + else if (ret == GLP_EBADB)
1.228 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.229 + xprintf("glp_simplex: initial basis is invalid\n");
1.230 + }
1.231 + else if (ret == GLP_ESING)
1.232 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.233 + xprintf("glp_simplex: initial basis is singular\n");
1.234 + }
1.235 + else if (ret == GLP_ECOND)
1.236 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.237 + xprintf(
1.238 + "glp_simplex: initial basis is ill-conditioned\n");
1.239 + }
1.240 + else
1.241 + xassert(ret != ret);
1.242 + if (ret != 0) goto done;
1.243 + }
1.244 + if (parm->meth == GLP_PRIMAL)
1.245 + ret = spx_primal(P, parm);
1.246 + else if (parm->meth == GLP_DUALP)
1.247 + { ret = spx_dual(P, parm);
1.248 + if (ret == GLP_EFAIL && P->valid)
1.249 + ret = spx_primal(P, parm);
1.250 + }
1.251 + else if (parm->meth == GLP_DUAL)
1.252 + ret = spx_dual(P, parm);
1.253 + else
1.254 + xassert(parm != parm);
1.255 +done: return ret;
1.256 +}
1.257 +
1.258 +static int preprocess_and_solve_lp(glp_prob *P, const glp_smcp *parm)
1.259 +{ /* solve LP using the preprocessor */
1.260 + NPP *npp;
1.261 + glp_prob *lp = NULL;
1.262 + glp_bfcp bfcp;
1.263 + int ret;
1.264 + if (parm->msg_lev >= GLP_MSG_ALL)
1.265 + xprintf("Preprocessing...\n");
1.266 + /* create preprocessor workspace */
1.267 + npp = npp_create_wksp();
1.268 + /* load original problem into the preprocessor workspace */
1.269 + npp_load_prob(npp, P, GLP_OFF, GLP_SOL, GLP_OFF);
1.270 + /* process LP prior to applying primal/dual simplex method */
1.271 + ret = npp_simplex(npp, parm);
1.272 + if (ret == 0)
1.273 + ;
1.274 + else if (ret == GLP_ENOPFS)
1.275 + { if (parm->msg_lev >= GLP_MSG_ALL)
1.276 + xprintf("PROBLEM HAS NO PRIMAL FEASIBLE SOLUTION\n");
1.277 + }
1.278 + else if (ret == GLP_ENODFS)
1.279 + { if (parm->msg_lev >= GLP_MSG_ALL)
1.280 + xprintf("PROBLEM HAS NO DUAL FEASIBLE SOLUTION\n");
1.281 + }
1.282 + else
1.283 + xassert(ret != ret);
1.284 + if (ret != 0) goto done;
1.285 + /* build transformed LP */
1.286 + lp = glp_create_prob();
1.287 + npp_build_prob(npp, lp);
1.288 + /* if the transformed LP is empty, it has empty solution, which
1.289 + is optimal */
1.290 + if (lp->m == 0 && lp->n == 0)
1.291 + { lp->pbs_stat = lp->dbs_stat = GLP_FEAS;
1.292 + lp->obj_val = lp->c0;
1.293 + if (parm->msg_lev >= GLP_MSG_ON && parm->out_dly == 0)
1.294 + { xprintf("~%6d: obj = %17.9e infeas = %10.3e\n", P->it_cnt,
1.295 + lp->obj_val, 0.0);
1.296 + }
1.297 + if (parm->msg_lev >= GLP_MSG_ALL)
1.298 + xprintf("OPTIMAL SOLUTION FOUND BY LP PREPROCESSOR\n");
1.299 + goto post;
1.300 + }
1.301 + if (parm->msg_lev >= GLP_MSG_ALL)
1.302 + { xprintf("%d row%s, %d column%s, %d non-zero%s\n",
1.303 + lp->m, lp->m == 1 ? "" : "s", lp->n, lp->n == 1 ? "" : "s",
1.304 + lp->nnz, lp->nnz == 1 ? "" : "s");
1.305 + }
1.306 + /* inherit basis factorization control parameters */
1.307 + glp_get_bfcp(P, &bfcp);
1.308 + glp_set_bfcp(lp, &bfcp);
1.309 + /* scale the transformed problem */
1.310 + { ENV *env = get_env_ptr();
1.311 + int term_out = env->term_out;
1.312 + if (!term_out || parm->msg_lev < GLP_MSG_ALL)
1.313 + env->term_out = GLP_OFF;
1.314 + else
1.315 + env->term_out = GLP_ON;
1.316 + glp_scale_prob(lp, GLP_SF_AUTO);
1.317 + env->term_out = term_out;
1.318 + }
1.319 + /* build advanced initial basis */
1.320 + { ENV *env = get_env_ptr();
1.321 + int term_out = env->term_out;
1.322 + if (!term_out || parm->msg_lev < GLP_MSG_ALL)
1.323 + env->term_out = GLP_OFF;
1.324 + else
1.325 + env->term_out = GLP_ON;
1.326 + glp_adv_basis(lp, 0);
1.327 + env->term_out = term_out;
1.328 + }
1.329 + /* solve the transformed LP */
1.330 + lp->it_cnt = P->it_cnt;
1.331 + ret = solve_lp(lp, parm);
1.332 + P->it_cnt = lp->it_cnt;
1.333 + /* only optimal solution can be postprocessed */
1.334 + if (!(ret == 0 && lp->pbs_stat == GLP_FEAS && lp->dbs_stat ==
1.335 + GLP_FEAS))
1.336 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.337 + xprintf("glp_simplex: unable to recover undefined or non-op"
1.338 + "timal solution\n");
1.339 + if (ret == 0)
1.340 + { if (lp->pbs_stat == GLP_NOFEAS)
1.341 + ret = GLP_ENOPFS;
1.342 + else if (lp->dbs_stat == GLP_NOFEAS)
1.343 + ret = GLP_ENODFS;
1.344 + else
1.345 + xassert(lp != lp);
1.346 + }
1.347 + goto done;
1.348 + }
1.349 +post: /* postprocess solution from the transformed LP */
1.350 + npp_postprocess(npp, lp);
1.351 + /* the transformed LP is no longer needed */
1.352 + glp_delete_prob(lp), lp = NULL;
1.353 + /* store solution to the original problem */
1.354 + npp_unload_sol(npp, P);
1.355 + /* the original LP has been successfully solved */
1.356 + ret = 0;
1.357 +done: /* delete the transformed LP, if it exists */
1.358 + if (lp != NULL) glp_delete_prob(lp);
1.359 + /* delete preprocessor workspace */
1.360 + npp_delete_wksp(npp);
1.361 + return ret;
1.362 +}
1.363 +
1.364 +int glp_simplex(glp_prob *P, const glp_smcp *parm)
1.365 +{ /* solve LP problem with the simplex method */
1.366 + glp_smcp _parm;
1.367 + int i, j, ret;
1.368 + /* check problem object */
1.369 + if (P == NULL || P->magic != GLP_PROB_MAGIC)
1.370 + xerror("glp_simplex: P = %p; invalid problem object\n", P);
1.371 + if (P->tree != NULL && P->tree->reason != 0)
1.372 + xerror("glp_simplex: operation not allowed\n");
1.373 + /* check control parameters */
1.374 + if (parm == NULL)
1.375 + parm = &_parm, glp_init_smcp((glp_smcp *)parm);
1.376 + if (!(parm->msg_lev == GLP_MSG_OFF ||
1.377 + parm->msg_lev == GLP_MSG_ERR ||
1.378 + parm->msg_lev == GLP_MSG_ON ||
1.379 + parm->msg_lev == GLP_MSG_ALL ||
1.380 + parm->msg_lev == GLP_MSG_DBG))
1.381 + xerror("glp_simplex: msg_lev = %d; invalid parameter\n",
1.382 + parm->msg_lev);
1.383 + if (!(parm->meth == GLP_PRIMAL ||
1.384 + parm->meth == GLP_DUALP ||
1.385 + parm->meth == GLP_DUAL))
1.386 + xerror("glp_simplex: meth = %d; invalid parameter\n",
1.387 + parm->meth);
1.388 + if (!(parm->pricing == GLP_PT_STD ||
1.389 + parm->pricing == GLP_PT_PSE))
1.390 + xerror("glp_simplex: pricing = %d; invalid parameter\n",
1.391 + parm->pricing);
1.392 + if (!(parm->r_test == GLP_RT_STD ||
1.393 + parm->r_test == GLP_RT_HAR))
1.394 + xerror("glp_simplex: r_test = %d; invalid parameter\n",
1.395 + parm->r_test);
1.396 + if (!(0.0 < parm->tol_bnd && parm->tol_bnd < 1.0))
1.397 + xerror("glp_simplex: tol_bnd = %g; invalid parameter\n",
1.398 + parm->tol_bnd);
1.399 + if (!(0.0 < parm->tol_dj && parm->tol_dj < 1.0))
1.400 + xerror("glp_simplex: tol_dj = %g; invalid parameter\n",
1.401 + parm->tol_dj);
1.402 + if (!(0.0 < parm->tol_piv && parm->tol_piv < 1.0))
1.403 + xerror("glp_simplex: tol_piv = %g; invalid parameter\n",
1.404 + parm->tol_piv);
1.405 + if (parm->it_lim < 0)
1.406 + xerror("glp_simplex: it_lim = %d; invalid parameter\n",
1.407 + parm->it_lim);
1.408 + if (parm->tm_lim < 0)
1.409 + xerror("glp_simplex: tm_lim = %d; invalid parameter\n",
1.410 + parm->tm_lim);
1.411 + if (parm->out_frq < 1)
1.412 + xerror("glp_simplex: out_frq = %d; invalid parameter\n",
1.413 + parm->out_frq);
1.414 + if (parm->out_dly < 0)
1.415 + xerror("glp_simplex: out_dly = %d; invalid parameter\n",
1.416 + parm->out_dly);
1.417 + if (!(parm->presolve == GLP_ON || parm->presolve == GLP_OFF))
1.418 + xerror("glp_simplex: presolve = %d; invalid parameter\n",
1.419 + parm->presolve);
1.420 + /* basic solution is currently undefined */
1.421 + P->pbs_stat = P->dbs_stat = GLP_UNDEF;
1.422 + P->obj_val = 0.0;
1.423 + P->some = 0;
1.424 + /* check bounds of double-bounded variables */
1.425 + for (i = 1; i <= P->m; i++)
1.426 + { GLPROW *row = P->row[i];
1.427 + if (row->type == GLP_DB && row->lb >= row->ub)
1.428 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.429 + xprintf("glp_simplex: row %d: lb = %g, ub = %g; incorrec"
1.430 + "t bounds\n", i, row->lb, row->ub);
1.431 + ret = GLP_EBOUND;
1.432 + goto done;
1.433 + }
1.434 + }
1.435 + for (j = 1; j <= P->n; j++)
1.436 + { GLPCOL *col = P->col[j];
1.437 + if (col->type == GLP_DB && col->lb >= col->ub)
1.438 + { if (parm->msg_lev >= GLP_MSG_ERR)
1.439 + xprintf("glp_simplex: column %d: lb = %g, ub = %g; incor"
1.440 + "rect bounds\n", j, col->lb, col->ub);
1.441 + ret = GLP_EBOUND;
1.442 + goto done;
1.443 + }
1.444 + }
1.445 + /* solve LP problem */
1.446 + if (parm->msg_lev >= GLP_MSG_ALL)
1.447 + { xprintf("GLPK Simplex Optimizer, v%s\n", glp_version());
1.448 + xprintf("%d row%s, %d column%s, %d non-zero%s\n",
1.449 + P->m, P->m == 1 ? "" : "s", P->n, P->n == 1 ? "" : "s",
1.450 + P->nnz, P->nnz == 1 ? "" : "s");
1.451 + }
1.452 + if (P->nnz == 0)
1.453 + trivial_lp(P, parm), ret = 0;
1.454 + else if (!parm->presolve)
1.455 + ret = solve_lp(P, parm);
1.456 + else
1.457 + ret = preprocess_and_solve_lp(P, parm);
1.458 +done: /* return to the application program */
1.459 + return ret;
1.460 +}
1.461 +
1.462 +/***********************************************************************
1.463 +* NAME
1.464 +*
1.465 +* glp_init_smcp - initialize simplex method control parameters
1.466 +*
1.467 +* SYNOPSIS
1.468 +*
1.469 +* void glp_init_smcp(glp_smcp *parm);
1.470 +*
1.471 +* DESCRIPTION
1.472 +*
1.473 +* The routine glp_init_smcp initializes control parameters, which are
1.474 +* used by the simplex solver, with default values.
1.475 +*
1.476 +* Default values of the control parameters are stored in a glp_smcp
1.477 +* structure, which the parameter parm points to. */
1.478 +
1.479 +void glp_init_smcp(glp_smcp *parm)
1.480 +{ parm->msg_lev = GLP_MSG_ALL;
1.481 + parm->meth = GLP_PRIMAL;
1.482 + parm->pricing = GLP_PT_PSE;
1.483 + parm->r_test = GLP_RT_HAR;
1.484 + parm->tol_bnd = 1e-7;
1.485 + parm->tol_dj = 1e-7;
1.486 + parm->tol_piv = 1e-10;
1.487 + parm->obj_ll = -DBL_MAX;
1.488 + parm->obj_ul = +DBL_MAX;
1.489 + parm->it_lim = INT_MAX;
1.490 + parm->tm_lim = INT_MAX;
1.491 + parm->out_frq = 500;
1.492 + parm->out_dly = 0;
1.493 + parm->presolve = GLP_OFF;
1.494 + return;
1.495 +}
1.496 +
1.497 +/***********************************************************************
1.498 +* NAME
1.499 +*
1.500 +* glp_get_status - retrieve generic status of basic solution
1.501 +*
1.502 +* SYNOPSIS
1.503 +*
1.504 +* int glp_get_status(glp_prob *lp);
1.505 +*
1.506 +* RETURNS
1.507 +*
1.508 +* The routine glp_get_status reports the generic status of the basic
1.509 +* solution for the specified problem object as follows:
1.510 +*
1.511 +* GLP_OPT - solution is optimal;
1.512 +* GLP_FEAS - solution is feasible;
1.513 +* GLP_INFEAS - solution is infeasible;
1.514 +* GLP_NOFEAS - problem has no feasible solution;
1.515 +* GLP_UNBND - problem has unbounded solution;
1.516 +* GLP_UNDEF - solution is undefined. */
1.517 +
1.518 +int glp_get_status(glp_prob *lp)
1.519 +{ int status;
1.520 + status = glp_get_prim_stat(lp);
1.521 + switch (status)
1.522 + { case GLP_FEAS:
1.523 + switch (glp_get_dual_stat(lp))
1.524 + { case GLP_FEAS:
1.525 + status = GLP_OPT;
1.526 + break;
1.527 + case GLP_NOFEAS:
1.528 + status = GLP_UNBND;
1.529 + break;
1.530 + case GLP_UNDEF:
1.531 + case GLP_INFEAS:
1.532 + status = status;
1.533 + break;
1.534 + default:
1.535 + xassert(lp != lp);
1.536 + }
1.537 + break;
1.538 + case GLP_UNDEF:
1.539 + case GLP_INFEAS:
1.540 + case GLP_NOFEAS:
1.541 + status = status;
1.542 + break;
1.543 + default:
1.544 + xassert(lp != lp);
1.545 + }
1.546 + return status;
1.547 +}
1.548 +
1.549 +/***********************************************************************
1.550 +* NAME
1.551 +*
1.552 +* glp_get_prim_stat - retrieve status of primal basic solution
1.553 +*
1.554 +* SYNOPSIS
1.555 +*
1.556 +* int glp_get_prim_stat(glp_prob *lp);
1.557 +*
1.558 +* RETURNS
1.559 +*
1.560 +* The routine glp_get_prim_stat reports the status of the primal basic
1.561 +* solution for the specified problem object as follows:
1.562 +*
1.563 +* GLP_UNDEF - primal solution is undefined;
1.564 +* GLP_FEAS - primal solution is feasible;
1.565 +* GLP_INFEAS - primal solution is infeasible;
1.566 +* GLP_NOFEAS - no primal feasible solution exists. */
1.567 +
1.568 +int glp_get_prim_stat(glp_prob *lp)
1.569 +{ int pbs_stat = lp->pbs_stat;
1.570 + return pbs_stat;
1.571 +}
1.572 +
1.573 +/***********************************************************************
1.574 +* NAME
1.575 +*
1.576 +* glp_get_dual_stat - retrieve status of dual basic solution
1.577 +*
1.578 +* SYNOPSIS
1.579 +*
1.580 +* int glp_get_dual_stat(glp_prob *lp);
1.581 +*
1.582 +* RETURNS
1.583 +*
1.584 +* The routine glp_get_dual_stat reports the status of the dual basic
1.585 +* solution for the specified problem object as follows:
1.586 +*
1.587 +* GLP_UNDEF - dual solution is undefined;
1.588 +* GLP_FEAS - dual solution is feasible;
1.589 +* GLP_INFEAS - dual solution is infeasible;
1.590 +* GLP_NOFEAS - no dual feasible solution exists. */
1.591 +
1.592 +int glp_get_dual_stat(glp_prob *lp)
1.593 +{ int dbs_stat = lp->dbs_stat;
1.594 + return dbs_stat;
1.595 +}
1.596 +
1.597 +/***********************************************************************
1.598 +* NAME
1.599 +*
1.600 +* glp_get_obj_val - retrieve objective value (basic solution)
1.601 +*
1.602 +* SYNOPSIS
1.603 +*
1.604 +* double glp_get_obj_val(glp_prob *lp);
1.605 +*
1.606 +* RETURNS
1.607 +*
1.608 +* The routine glp_get_obj_val returns value of the objective function
1.609 +* for basic solution. */
1.610 +
1.611 +double glp_get_obj_val(glp_prob *lp)
1.612 +{ /*struct LPXCPS *cps = lp->cps;*/
1.613 + double z;
1.614 + z = lp->obj_val;
1.615 + /*if (cps->round && fabs(z) < 1e-9) z = 0.0;*/
1.616 + return z;
1.617 +}
1.618 +
1.619 +/***********************************************************************
1.620 +* NAME
1.621 +*
1.622 +* glp_get_row_stat - retrieve row status
1.623 +*
1.624 +* SYNOPSIS
1.625 +*
1.626 +* int glp_get_row_stat(glp_prob *lp, int i);
1.627 +*
1.628 +* RETURNS
1.629 +*
1.630 +* The routine glp_get_row_stat returns current status assigned to the
1.631 +* auxiliary variable associated with i-th row as follows:
1.632 +*
1.633 +* GLP_BS - basic variable;
1.634 +* GLP_NL - non-basic variable on its lower bound;
1.635 +* GLP_NU - non-basic variable on its upper bound;
1.636 +* GLP_NF - non-basic free (unbounded) variable;
1.637 +* GLP_NS - non-basic fixed variable. */
1.638 +
1.639 +int glp_get_row_stat(glp_prob *lp, int i)
1.640 +{ if (!(1 <= i && i <= lp->m))
1.641 + xerror("glp_get_row_stat: i = %d; row number out of range\n",
1.642 + i);
1.643 + return lp->row[i]->stat;
1.644 +}
1.645 +
1.646 +/***********************************************************************
1.647 +* NAME
1.648 +*
1.649 +* glp_get_row_prim - retrieve row primal value (basic solution)
1.650 +*
1.651 +* SYNOPSIS
1.652 +*
1.653 +* double glp_get_row_prim(glp_prob *lp, int i);
1.654 +*
1.655 +* RETURNS
1.656 +*
1.657 +* The routine glp_get_row_prim returns primal value of the auxiliary
1.658 +* variable associated with i-th row. */
1.659 +
1.660 +double glp_get_row_prim(glp_prob *lp, int i)
1.661 +{ /*struct LPXCPS *cps = lp->cps;*/
1.662 + double prim;
1.663 + if (!(1 <= i && i <= lp->m))
1.664 + xerror("glp_get_row_prim: i = %d; row number out of range\n",
1.665 + i);
1.666 + prim = lp->row[i]->prim;
1.667 + /*if (cps->round && fabs(prim) < 1e-9) prim = 0.0;*/
1.668 + return prim;
1.669 +}
1.670 +
1.671 +/***********************************************************************
1.672 +* NAME
1.673 +*
1.674 +* glp_get_row_dual - retrieve row dual value (basic solution)
1.675 +*
1.676 +* SYNOPSIS
1.677 +*
1.678 +* double glp_get_row_dual(glp_prob *lp, int i);
1.679 +*
1.680 +* RETURNS
1.681 +*
1.682 +* The routine glp_get_row_dual returns dual value (i.e. reduced cost)
1.683 +* of the auxiliary variable associated with i-th row. */
1.684 +
1.685 +double glp_get_row_dual(glp_prob *lp, int i)
1.686 +{ /*struct LPXCPS *cps = lp->cps;*/
1.687 + double dual;
1.688 + if (!(1 <= i && i <= lp->m))
1.689 + xerror("glp_get_row_dual: i = %d; row number out of range\n",
1.690 + i);
1.691 + dual = lp->row[i]->dual;
1.692 + /*if (cps->round && fabs(dual) < 1e-9) dual = 0.0;*/
1.693 + return dual;
1.694 +}
1.695 +
1.696 +/***********************************************************************
1.697 +* NAME
1.698 +*
1.699 +* glp_get_col_stat - retrieve column status
1.700 +*
1.701 +* SYNOPSIS
1.702 +*
1.703 +* int glp_get_col_stat(glp_prob *lp, int j);
1.704 +*
1.705 +* RETURNS
1.706 +*
1.707 +* The routine glp_get_col_stat returns current status assigned to the
1.708 +* structural variable associated with j-th column as follows:
1.709 +*
1.710 +* GLP_BS - basic variable;
1.711 +* GLP_NL - non-basic variable on its lower bound;
1.712 +* GLP_NU - non-basic variable on its upper bound;
1.713 +* GLP_NF - non-basic free (unbounded) variable;
1.714 +* GLP_NS - non-basic fixed variable. */
1.715 +
1.716 +int glp_get_col_stat(glp_prob *lp, int j)
1.717 +{ if (!(1 <= j && j <= lp->n))
1.718 + xerror("glp_get_col_stat: j = %d; column number out of range\n"
1.719 + , j);
1.720 + return lp->col[j]->stat;
1.721 +}
1.722 +
1.723 +/***********************************************************************
1.724 +* NAME
1.725 +*
1.726 +* glp_get_col_prim - retrieve column primal value (basic solution)
1.727 +*
1.728 +* SYNOPSIS
1.729 +*
1.730 +* double glp_get_col_prim(glp_prob *lp, int j);
1.731 +*
1.732 +* RETURNS
1.733 +*
1.734 +* The routine glp_get_col_prim returns primal value of the structural
1.735 +* variable associated with j-th column. */
1.736 +
1.737 +double glp_get_col_prim(glp_prob *lp, int j)
1.738 +{ /*struct LPXCPS *cps = lp->cps;*/
1.739 + double prim;
1.740 + if (!(1 <= j && j <= lp->n))
1.741 + xerror("glp_get_col_prim: j = %d; column number out of range\n"
1.742 + , j);
1.743 + prim = lp->col[j]->prim;
1.744 + /*if (cps->round && fabs(prim) < 1e-9) prim = 0.0;*/
1.745 + return prim;
1.746 +}
1.747 +
1.748 +/***********************************************************************
1.749 +* NAME
1.750 +*
1.751 +* glp_get_col_dual - retrieve column dual value (basic solution)
1.752 +*
1.753 +* SYNOPSIS
1.754 +*
1.755 +* double glp_get_col_dual(glp_prob *lp, int j);
1.756 +*
1.757 +* RETURNS
1.758 +*
1.759 +* The routine glp_get_col_dual returns dual value (i.e. reduced cost)
1.760 +* of the structural variable associated with j-th column. */
1.761 +
1.762 +double glp_get_col_dual(glp_prob *lp, int j)
1.763 +{ /*struct LPXCPS *cps = lp->cps;*/
1.764 + double dual;
1.765 + if (!(1 <= j && j <= lp->n))
1.766 + xerror("glp_get_col_dual: j = %d; column number out of range\n"
1.767 + , j);
1.768 + dual = lp->col[j]->dual;
1.769 + /*if (cps->round && fabs(dual) < 1e-9) dual = 0.0;*/
1.770 + return dual;
1.771 +}
1.772 +
1.773 +/***********************************************************************
1.774 +* NAME
1.775 +*
1.776 +* glp_get_unbnd_ray - determine variable causing unboundedness
1.777 +*
1.778 +* SYNOPSIS
1.779 +*
1.780 +* int glp_get_unbnd_ray(glp_prob *lp);
1.781 +*
1.782 +* RETURNS
1.783 +*
1.784 +* The routine glp_get_unbnd_ray returns the number k of a variable,
1.785 +* which causes primal or dual unboundedness. If 1 <= k <= m, it is
1.786 +* k-th auxiliary variable, and if m+1 <= k <= m+n, it is (k-m)-th
1.787 +* structural variable, where m is the number of rows, n is the number
1.788 +* of columns in the problem object. If such variable is not defined,
1.789 +* the routine returns 0.
1.790 +*
1.791 +* COMMENTS
1.792 +*
1.793 +* If it is not exactly known which version of the simplex solver
1.794 +* detected unboundedness, i.e. whether the unboundedness is primal or
1.795 +* dual, it is sufficient to check the status of the variable reported
1.796 +* with the routine glp_get_row_stat or glp_get_col_stat. If the
1.797 +* variable is non-basic, the unboundedness is primal, otherwise, if
1.798 +* the variable is basic, the unboundedness is dual (the latter case
1.799 +* means that the problem has no primal feasible dolution). */
1.800 +
1.801 +int glp_get_unbnd_ray(glp_prob *lp)
1.802 +{ int k;
1.803 + k = lp->some;
1.804 + xassert(k >= 0);
1.805 + if (k > lp->m + lp->n) k = 0;
1.806 + return k;
1.807 +}
1.808 +
1.809 +/* eof */