/* glpapi06.c (simplex method routines) */ /*********************************************************************** * This code is part of GLPK (GNU Linear Programming Kit). * * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, * 2009, 2010 Andrew Makhorin, Department for Applied Informatics, * Moscow Aviation Institute, Moscow, Russia. All rights reserved. * E-mail: . * * GLPK is free software: you can redistribute it and/or modify it * under the terms of the GNU General Public License as published by * the Free Software Foundation, either version 3 of the License, or * (at your option) any later version. * * GLPK is distributed in the hope that it will be useful, but WITHOUT * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public * License for more details. * * You should have received a copy of the GNU General Public License * along with GLPK. If not, see . ***********************************************************************/ #include "glpios.h" #include "glpnpp.h" #include "glpspx.h" /*********************************************************************** * NAME * * glp_simplex - solve LP problem with the simplex method * * SYNOPSIS * * int glp_simplex(glp_prob *P, const glp_smcp *parm); * * DESCRIPTION * * The routine glp_simplex is a driver to the LP solver based on the * simplex method. This routine retrieves problem data from the * specified problem object, calls the solver to solve the problem * instance, and stores results of computations back into the problem * object. * * The simplex solver has a set of control parameters. Values of the * control parameters can be passed in a structure glp_smcp, which the * parameter parm points to. * * The parameter parm can be specified as NULL, in which case the LP * solver uses default settings. * * RETURNS * * 0 The LP problem instance has been successfully solved. This code * does not necessarily mean that the solver has found optimal * solution. It only means that the solution process was successful. * * GLP_EBADB * Unable to start the search, because the initial basis specified * in the problem object is invalid--the number of basic (auxiliary * and structural) variables is not the same as the number of rows in * the problem object. * * GLP_ESING * Unable to start the search, because the basis matrix correspodning * to the initial basis is singular within the working precision. * * GLP_ECOND * Unable to start the search, because the basis matrix correspodning * to the initial basis is ill-conditioned, i.e. its condition number * is too large. * * GLP_EBOUND * Unable to start the search, because some double-bounded variables * have incorrect bounds. * * GLP_EFAIL * The search was prematurely terminated due to the solver failure. * * GLP_EOBJLL * The search was prematurely terminated, because the objective * function being maximized has reached its lower limit and continues * decreasing (dual simplex only). * * GLP_EOBJUL * The search was prematurely terminated, because the objective * function being minimized has reached its upper limit and continues * increasing (dual simplex only). * * GLP_EITLIM * The search was prematurely terminated, because the simplex * iteration limit has been exceeded. * * GLP_ETMLIM * The search was prematurely terminated, because the time limit has * been exceeded. * * GLP_ENOPFS * The LP problem instance has no primal feasible solution (only if * the LP presolver is used). * * GLP_ENODFS * The LP problem instance has no dual feasible solution (only if the * LP presolver is used). */ static void trivial_lp(glp_prob *P, const glp_smcp *parm) { /* solve trivial LP which has empty constraint matrix */ GLPROW *row; GLPCOL *col; int i, j; double p_infeas, d_infeas, zeta; P->valid = 0; P->pbs_stat = P->dbs_stat = GLP_FEAS; P->obj_val = P->c0; P->some = 0; p_infeas = d_infeas = 0.0; /* make all auxiliary variables basic */ for (i = 1; i <= P->m; i++) { row = P->row[i]; row->stat = GLP_BS; row->prim = row->dual = 0.0; /* check primal feasibility */ if (row->type == GLP_LO || row->type == GLP_DB || row->type == GLP_FX) { /* row has lower bound */ if (row->lb > + parm->tol_bnd) { P->pbs_stat = GLP_NOFEAS; if (P->some == 0 && parm->meth != GLP_PRIMAL) P->some = i; } if (p_infeas < + row->lb) p_infeas = + row->lb; } if (row->type == GLP_UP || row->type == GLP_DB || row->type == GLP_FX) { /* row has upper bound */ if (row->ub < - parm->tol_bnd) { P->pbs_stat = GLP_NOFEAS; if (P->some == 0 && parm->meth != GLP_PRIMAL) P->some = i; } if (p_infeas < - row->ub) p_infeas = - row->ub; } } /* determine scale factor for the objective row */ zeta = 1.0; for (j = 1; j <= P->n; j++) { col = P->col[j]; if (zeta < fabs(col->coef)) zeta = fabs(col->coef); } zeta = (P->dir == GLP_MIN ? +1.0 : -1.0) / zeta; /* make all structural variables non-basic */ for (j = 1; j <= P->n; j++) { col = P->col[j]; if (col->type == GLP_FR) col->stat = GLP_NF, col->prim = 0.0; else if (col->type == GLP_LO) lo: col->stat = GLP_NL, col->prim = col->lb; else if (col->type == GLP_UP) up: col->stat = GLP_NU, col->prim = col->ub; else if (col->type == GLP_DB) { if (zeta * col->coef > 0.0) goto lo; else if (zeta * col->coef < 0.0) goto up; else if (fabs(col->lb) <= fabs(col->ub)) goto lo; else goto up; } else if (col->type == GLP_FX) col->stat = GLP_NS, col->prim = col->lb; col->dual = col->coef; P->obj_val += col->coef * col->prim; /* check dual feasibility */ if (col->type == GLP_FR || col->type == GLP_LO) { /* column has no upper bound */ if (zeta * col->dual < - parm->tol_dj) { P->dbs_stat = GLP_NOFEAS; if (P->some == 0 && parm->meth == GLP_PRIMAL) P->some = P->m + j; } if (d_infeas < - zeta * col->dual) d_infeas = - zeta * col->dual; } if (col->type == GLP_FR || col->type == GLP_UP) { /* column has no lower bound */ if (zeta * col->dual > + parm->tol_dj) { P->dbs_stat = GLP_NOFEAS; if (P->some == 0 && parm->meth == GLP_PRIMAL) P->some = P->m + j; } if (d_infeas < + zeta * col->dual) d_infeas = + zeta * col->dual; } } /* simulate the simplex solver output */ if (parm->msg_lev >= GLP_MSG_ON && parm->out_dly == 0) { xprintf("~%6d: obj = %17.9e infeas = %10.3e\n", P->it_cnt, P->obj_val, parm->meth == GLP_PRIMAL ? p_infeas : d_infeas); } if (parm->msg_lev >= GLP_MSG_ALL && parm->out_dly == 0) { if (P->pbs_stat == GLP_FEAS && P->dbs_stat == GLP_FEAS) xprintf("OPTIMAL SOLUTION FOUND\n"); else if (P->pbs_stat == GLP_NOFEAS) xprintf("PROBLEM HAS NO FEASIBLE SOLUTION\n"); else if (parm->meth == GLP_PRIMAL) xprintf("PROBLEM HAS UNBOUNDED SOLUTION\n"); else xprintf("PROBLEM HAS NO DUAL FEASIBLE SOLUTION\n"); } return; } static int solve_lp(glp_prob *P, const glp_smcp *parm) { /* solve LP directly without using the preprocessor */ int ret; if (!glp_bf_exists(P)) { ret = glp_factorize(P); if (ret == 0) ; else if (ret == GLP_EBADB) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf("glp_simplex: initial basis is invalid\n"); } else if (ret == GLP_ESING) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf("glp_simplex: initial basis is singular\n"); } else if (ret == GLP_ECOND) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf( "glp_simplex: initial basis is ill-conditioned\n"); } else xassert(ret != ret); if (ret != 0) goto done; } if (parm->meth == GLP_PRIMAL) ret = spx_primal(P, parm); else if (parm->meth == GLP_DUALP) { ret = spx_dual(P, parm); if (ret == GLP_EFAIL && P->valid) ret = spx_primal(P, parm); } else if (parm->meth == GLP_DUAL) ret = spx_dual(P, parm); else xassert(parm != parm); done: return ret; } static int preprocess_and_solve_lp(glp_prob *P, const glp_smcp *parm) { /* solve LP using the preprocessor */ NPP *npp; glp_prob *lp = NULL; glp_bfcp bfcp; int ret; if (parm->msg_lev >= GLP_MSG_ALL) xprintf("Preprocessing...\n"); /* create preprocessor workspace */ npp = npp_create_wksp(); /* load original problem into the preprocessor workspace */ npp_load_prob(npp, P, GLP_OFF, GLP_SOL, GLP_OFF); /* process LP prior to applying primal/dual simplex method */ ret = npp_simplex(npp, parm); if (ret == 0) ; else if (ret == GLP_ENOPFS) { if (parm->msg_lev >= GLP_MSG_ALL) xprintf("PROBLEM HAS NO PRIMAL FEASIBLE SOLUTION\n"); } else if (ret == GLP_ENODFS) { if (parm->msg_lev >= GLP_MSG_ALL) xprintf("PROBLEM HAS NO DUAL FEASIBLE SOLUTION\n"); } else xassert(ret != ret); if (ret != 0) goto done; /* build transformed LP */ lp = glp_create_prob(); npp_build_prob(npp, lp); /* if the transformed LP is empty, it has empty solution, which is optimal */ if (lp->m == 0 && lp->n == 0) { lp->pbs_stat = lp->dbs_stat = GLP_FEAS; lp->obj_val = lp->c0; if (parm->msg_lev >= GLP_MSG_ON && parm->out_dly == 0) { xprintf("~%6d: obj = %17.9e infeas = %10.3e\n", P->it_cnt, lp->obj_val, 0.0); } if (parm->msg_lev >= GLP_MSG_ALL) xprintf("OPTIMAL SOLUTION FOUND BY LP PREPROCESSOR\n"); goto post; } if (parm->msg_lev >= GLP_MSG_ALL) { xprintf("%d row%s, %d column%s, %d non-zero%s\n", lp->m, lp->m == 1 ? "" : "s", lp->n, lp->n == 1 ? "" : "s", lp->nnz, lp->nnz == 1 ? "" : "s"); } /* inherit basis factorization control parameters */ glp_get_bfcp(P, &bfcp); glp_set_bfcp(lp, &bfcp); /* scale the transformed problem */ { ENV *env = get_env_ptr(); int term_out = env->term_out; if (!term_out || parm->msg_lev < GLP_MSG_ALL) env->term_out = GLP_OFF; else env->term_out = GLP_ON; glp_scale_prob(lp, GLP_SF_AUTO); env->term_out = term_out; } /* build advanced initial basis */ { ENV *env = get_env_ptr(); int term_out = env->term_out; if (!term_out || parm->msg_lev < GLP_MSG_ALL) env->term_out = GLP_OFF; else env->term_out = GLP_ON; glp_adv_basis(lp, 0); env->term_out = term_out; } /* solve the transformed LP */ lp->it_cnt = P->it_cnt; ret = solve_lp(lp, parm); P->it_cnt = lp->it_cnt; /* only optimal solution can be postprocessed */ if (!(ret == 0 && lp->pbs_stat == GLP_FEAS && lp->dbs_stat == GLP_FEAS)) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf("glp_simplex: unable to recover undefined or non-op" "timal solution\n"); if (ret == 0) { if (lp->pbs_stat == GLP_NOFEAS) ret = GLP_ENOPFS; else if (lp->dbs_stat == GLP_NOFEAS) ret = GLP_ENODFS; else xassert(lp != lp); } goto done; } post: /* postprocess solution from the transformed LP */ npp_postprocess(npp, lp); /* the transformed LP is no longer needed */ glp_delete_prob(lp), lp = NULL; /* store solution to the original problem */ npp_unload_sol(npp, P); /* the original LP has been successfully solved */ ret = 0; done: /* delete the transformed LP, if it exists */ if (lp != NULL) glp_delete_prob(lp); /* delete preprocessor workspace */ npp_delete_wksp(npp); return ret; } int glp_simplex(glp_prob *P, const glp_smcp *parm) { /* solve LP problem with the simplex method */ glp_smcp _parm; int i, j, ret; /* check problem object */ if (P == NULL || P->magic != GLP_PROB_MAGIC) xerror("glp_simplex: P = %p; invalid problem object\n", P); if (P->tree != NULL && P->tree->reason != 0) xerror("glp_simplex: operation not allowed\n"); /* check control parameters */ if (parm == NULL) parm = &_parm, glp_init_smcp((glp_smcp *)parm); if (!(parm->msg_lev == GLP_MSG_OFF || parm->msg_lev == GLP_MSG_ERR || parm->msg_lev == GLP_MSG_ON || parm->msg_lev == GLP_MSG_ALL || parm->msg_lev == GLP_MSG_DBG)) xerror("glp_simplex: msg_lev = %d; invalid parameter\n", parm->msg_lev); if (!(parm->meth == GLP_PRIMAL || parm->meth == GLP_DUALP || parm->meth == GLP_DUAL)) xerror("glp_simplex: meth = %d; invalid parameter\n", parm->meth); if (!(parm->pricing == GLP_PT_STD || parm->pricing == GLP_PT_PSE)) xerror("glp_simplex: pricing = %d; invalid parameter\n", parm->pricing); if (!(parm->r_test == GLP_RT_STD || parm->r_test == GLP_RT_HAR)) xerror("glp_simplex: r_test = %d; invalid parameter\n", parm->r_test); if (!(0.0 < parm->tol_bnd && parm->tol_bnd < 1.0)) xerror("glp_simplex: tol_bnd = %g; invalid parameter\n", parm->tol_bnd); if (!(0.0 < parm->tol_dj && parm->tol_dj < 1.0)) xerror("glp_simplex: tol_dj = %g; invalid parameter\n", parm->tol_dj); if (!(0.0 < parm->tol_piv && parm->tol_piv < 1.0)) xerror("glp_simplex: tol_piv = %g; invalid parameter\n", parm->tol_piv); if (parm->it_lim < 0) xerror("glp_simplex: it_lim = %d; invalid parameter\n", parm->it_lim); if (parm->tm_lim < 0) xerror("glp_simplex: tm_lim = %d; invalid parameter\n", parm->tm_lim); if (parm->out_frq < 1) xerror("glp_simplex: out_frq = %d; invalid parameter\n", parm->out_frq); if (parm->out_dly < 0) xerror("glp_simplex: out_dly = %d; invalid parameter\n", parm->out_dly); if (!(parm->presolve == GLP_ON || parm->presolve == GLP_OFF)) xerror("glp_simplex: presolve = %d; invalid parameter\n", parm->presolve); /* basic solution is currently undefined */ P->pbs_stat = P->dbs_stat = GLP_UNDEF; P->obj_val = 0.0; P->some = 0; /* check bounds of double-bounded variables */ for (i = 1; i <= P->m; i++) { GLPROW *row = P->row[i]; if (row->type == GLP_DB && row->lb >= row->ub) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf("glp_simplex: row %d: lb = %g, ub = %g; incorrec" "t bounds\n", i, row->lb, row->ub); ret = GLP_EBOUND; goto done; } } for (j = 1; j <= P->n; j++) { GLPCOL *col = P->col[j]; if (col->type == GLP_DB && col->lb >= col->ub) { if (parm->msg_lev >= GLP_MSG_ERR) xprintf("glp_simplex: column %d: lb = %g, ub = %g; incor" "rect bounds\n", j, col->lb, col->ub); ret = GLP_EBOUND; goto done; } } /* solve LP problem */ if (parm->msg_lev >= GLP_MSG_ALL) { xprintf("GLPK Simplex Optimizer, v%s\n", glp_version()); xprintf("%d row%s, %d column%s, %d non-zero%s\n", P->m, P->m == 1 ? "" : "s", P->n, P->n == 1 ? "" : "s", P->nnz, P->nnz == 1 ? "" : "s"); } if (P->nnz == 0) trivial_lp(P, parm), ret = 0; else if (!parm->presolve) ret = solve_lp(P, parm); else ret = preprocess_and_solve_lp(P, parm); done: /* return to the application program */ return ret; } /*********************************************************************** * NAME * * glp_init_smcp - initialize simplex method control parameters * * SYNOPSIS * * void glp_init_smcp(glp_smcp *parm); * * DESCRIPTION * * The routine glp_init_smcp initializes control parameters, which are * used by the simplex solver, with default values. * * Default values of the control parameters are stored in a glp_smcp * structure, which the parameter parm points to. */ void glp_init_smcp(glp_smcp *parm) { parm->msg_lev = GLP_MSG_ALL; parm->meth = GLP_PRIMAL; parm->pricing = GLP_PT_PSE; parm->r_test = GLP_RT_HAR; parm->tol_bnd = 1e-7; parm->tol_dj = 1e-7; parm->tol_piv = 1e-10; parm->obj_ll = -DBL_MAX; parm->obj_ul = +DBL_MAX; parm->it_lim = INT_MAX; parm->tm_lim = INT_MAX; parm->out_frq = 500; parm->out_dly = 0; parm->presolve = GLP_OFF; return; } /*********************************************************************** * NAME * * glp_get_status - retrieve generic status of basic solution * * SYNOPSIS * * int glp_get_status(glp_prob *lp); * * RETURNS * * The routine glp_get_status reports the generic status of the basic * solution for the specified problem object as follows: * * GLP_OPT - solution is optimal; * GLP_FEAS - solution is feasible; * GLP_INFEAS - solution is infeasible; * GLP_NOFEAS - problem has no feasible solution; * GLP_UNBND - problem has unbounded solution; * GLP_UNDEF - solution is undefined. */ int glp_get_status(glp_prob *lp) { int status; status = glp_get_prim_stat(lp); switch (status) { case GLP_FEAS: switch (glp_get_dual_stat(lp)) { case GLP_FEAS: status = GLP_OPT; break; case GLP_NOFEAS: status = GLP_UNBND; break; case GLP_UNDEF: case GLP_INFEAS: status = status; break; default: xassert(lp != lp); } break; case GLP_UNDEF: case GLP_INFEAS: case GLP_NOFEAS: status = status; break; default: xassert(lp != lp); } return status; } /*********************************************************************** * NAME * * glp_get_prim_stat - retrieve status of primal basic solution * * SYNOPSIS * * int glp_get_prim_stat(glp_prob *lp); * * RETURNS * * The routine glp_get_prim_stat reports the status of the primal basic * solution for the specified problem object as follows: * * GLP_UNDEF - primal solution is undefined; * GLP_FEAS - primal solution is feasible; * GLP_INFEAS - primal solution is infeasible; * GLP_NOFEAS - no primal feasible solution exists. */ int glp_get_prim_stat(glp_prob *lp) { int pbs_stat = lp->pbs_stat; return pbs_stat; } /*********************************************************************** * NAME * * glp_get_dual_stat - retrieve status of dual basic solution * * SYNOPSIS * * int glp_get_dual_stat(glp_prob *lp); * * RETURNS * * The routine glp_get_dual_stat reports the status of the dual basic * solution for the specified problem object as follows: * * GLP_UNDEF - dual solution is undefined; * GLP_FEAS - dual solution is feasible; * GLP_INFEAS - dual solution is infeasible; * GLP_NOFEAS - no dual feasible solution exists. */ int glp_get_dual_stat(glp_prob *lp) { int dbs_stat = lp->dbs_stat; return dbs_stat; } /*********************************************************************** * NAME * * glp_get_obj_val - retrieve objective value (basic solution) * * SYNOPSIS * * double glp_get_obj_val(glp_prob *lp); * * RETURNS * * The routine glp_get_obj_val returns value of the objective function * for basic solution. */ double glp_get_obj_val(glp_prob *lp) { /*struct LPXCPS *cps = lp->cps;*/ double z; z = lp->obj_val; /*if (cps->round && fabs(z) < 1e-9) z = 0.0;*/ return z; } /*********************************************************************** * NAME * * glp_get_row_stat - retrieve row status * * SYNOPSIS * * int glp_get_row_stat(glp_prob *lp, int i); * * RETURNS * * The routine glp_get_row_stat returns current status assigned to the * auxiliary variable associated with i-th row as follows: * * GLP_BS - basic variable; * GLP_NL - non-basic variable on its lower bound; * GLP_NU - non-basic variable on its upper bound; * GLP_NF - non-basic free (unbounded) variable; * GLP_NS - non-basic fixed variable. */ int glp_get_row_stat(glp_prob *lp, int i) { if (!(1 <= i && i <= lp->m)) xerror("glp_get_row_stat: i = %d; row number out of range\n", i); return lp->row[i]->stat; } /*********************************************************************** * NAME * * glp_get_row_prim - retrieve row primal value (basic solution) * * SYNOPSIS * * double glp_get_row_prim(glp_prob *lp, int i); * * RETURNS * * The routine glp_get_row_prim returns primal value of the auxiliary * variable associated with i-th row. */ double glp_get_row_prim(glp_prob *lp, int i) { /*struct LPXCPS *cps = lp->cps;*/ double prim; if (!(1 <= i && i <= lp->m)) xerror("glp_get_row_prim: i = %d; row number out of range\n", i); prim = lp->row[i]->prim; /*if (cps->round && fabs(prim) < 1e-9) prim = 0.0;*/ return prim; } /*********************************************************************** * NAME * * glp_get_row_dual - retrieve row dual value (basic solution) * * SYNOPSIS * * double glp_get_row_dual(glp_prob *lp, int i); * * RETURNS * * The routine glp_get_row_dual returns dual value (i.e. reduced cost) * of the auxiliary variable associated with i-th row. */ double glp_get_row_dual(glp_prob *lp, int i) { /*struct LPXCPS *cps = lp->cps;*/ double dual; if (!(1 <= i && i <= lp->m)) xerror("glp_get_row_dual: i = %d; row number out of range\n", i); dual = lp->row[i]->dual; /*if (cps->round && fabs(dual) < 1e-9) dual = 0.0;*/ return dual; } /*********************************************************************** * NAME * * glp_get_col_stat - retrieve column status * * SYNOPSIS * * int glp_get_col_stat(glp_prob *lp, int j); * * RETURNS * * The routine glp_get_col_stat returns current status assigned to the * structural variable associated with j-th column as follows: * * GLP_BS - basic variable; * GLP_NL - non-basic variable on its lower bound; * GLP_NU - non-basic variable on its upper bound; * GLP_NF - non-basic free (unbounded) variable; * GLP_NS - non-basic fixed variable. */ int glp_get_col_stat(glp_prob *lp, int j) { if (!(1 <= j && j <= lp->n)) xerror("glp_get_col_stat: j = %d; column number out of range\n" , j); return lp->col[j]->stat; } /*********************************************************************** * NAME * * glp_get_col_prim - retrieve column primal value (basic solution) * * SYNOPSIS * * double glp_get_col_prim(glp_prob *lp, int j); * * RETURNS * * The routine glp_get_col_prim returns primal value of the structural * variable associated with j-th column. */ double glp_get_col_prim(glp_prob *lp, int j) { /*struct LPXCPS *cps = lp->cps;*/ double prim; if (!(1 <= j && j <= lp->n)) xerror("glp_get_col_prim: j = %d; column number out of range\n" , j); prim = lp->col[j]->prim; /*if (cps->round && fabs(prim) < 1e-9) prim = 0.0;*/ return prim; } /*********************************************************************** * NAME * * glp_get_col_dual - retrieve column dual value (basic solution) * * SYNOPSIS * * double glp_get_col_dual(glp_prob *lp, int j); * * RETURNS * * The routine glp_get_col_dual returns dual value (i.e. reduced cost) * of the structural variable associated with j-th column. */ double glp_get_col_dual(glp_prob *lp, int j) { /*struct LPXCPS *cps = lp->cps;*/ double dual; if (!(1 <= j && j <= lp->n)) xerror("glp_get_col_dual: j = %d; column number out of range\n" , j); dual = lp->col[j]->dual; /*if (cps->round && fabs(dual) < 1e-9) dual = 0.0;*/ return dual; } /*********************************************************************** * NAME * * glp_get_unbnd_ray - determine variable causing unboundedness * * SYNOPSIS * * int glp_get_unbnd_ray(glp_prob *lp); * * RETURNS * * The routine glp_get_unbnd_ray returns the number k of a variable, * which causes primal or dual unboundedness. If 1 <= k <= m, it is * k-th auxiliary variable, and if m+1 <= k <= m+n, it is (k-m)-th * structural variable, where m is the number of rows, n is the number * of columns in the problem object. If such variable is not defined, * the routine returns 0. * * COMMENTS * * If it is not exactly known which version of the simplex solver * detected unboundedness, i.e. whether the unboundedness is primal or * dual, it is sufficient to check the status of the variable reported * with the routine glp_get_row_stat or glp_get_col_stat. If the * variable is non-basic, the unboundedness is primal, otherwise, if * the variable is basic, the unboundedness is dual (the latter case * means that the problem has no primal feasible dolution). */ int glp_get_unbnd_ray(glp_prob *lp) { int k; k = lp->some; xassert(k >= 0); if (k > lp->m + lp->n) k = 0; return k; } /* eof */