glplpx06.c

著名的大规模线性规划求解器源码GLPK.C语言版本,可以修剪.内有详细帮助文档.

      xfree(a);
      return len;
}

/*----------------------------------------------------------------------
-- lpx_prim_ratio_test - perform primal ratio test.
--
-- *Synopsis*
--
-- #include "glplpx.h"
-- int lpx_prim_ratio_test(LPX *lp, int len, int ind[], double val[],
--    int how, double tol);
--
-- *Description*
--
-- The routine lpx_prim_ratio_test performs the primal ratio test for
-- an explicitly specified column of the simplex table.
--
-- The primal basic solution associated with an LP problem object,
-- which the parameter lp points to, should be feasible. No components
-- of the LP problem object are changed by the routine.
--
-- The explicitly specified column of the simplex table shows how the
-- basic variables xB depend on some non-basic variable y (which is not
-- necessarily presented in the problem object):
--
--    xB[1] = ... + alfa[1]*y + ...
--    xB[2] = ... + alfa[2]*y + ...                                  (*)
--       .  .  .  .  .  .  .  .
--    xB[m] = ... + alfa[m]*y + ...
--
-- The column (*) is specifed on entry to the routine using the sparse
-- format. Ordinal numbers of basic variables xB[i] should be placed in
-- locations ind[1], ..., ind[len], where ordinal number 1 to m denote
-- auxiliary variables, and ordinal numbers m+1 to m+n denote structural
-- variables. The corresponding non-zero coefficients alfa[i] should be
-- placed in locations val[1], ..., val[len]. The arrays ind and val are
-- not changed on exit.
--
-- The parameter how specifies in which direction the variable y changes
-- on entering the basis: +1 means increasing, -1 means decreasing.
--
-- The parameter tol is a relative tolerance (small positive number)
-- used by the routine to skip small alfa[i] of the column (*).
--
-- The routine determines the ordinal number of some basic variable
-- (specified in ind[1], ..., ind[len]), which should leave the basis
-- instead the variable y in order to keep primal feasibility, and
-- returns it on exit. If the choice cannot be made (i.e. if the
-- adjacent basic solution is primal unbounded), the routine returns
-- zero.
--
-- *Note*
--
-- If the non-basic variable y is presented in the LP problem object,
-- the column (*) can be computed using the routine lpx_eval_tab_col.
-- Otherwise it can be computed using the routine lpx_transform_col.
--
-- *Returns*
--
-- The routine lpx_prim_ratio_test returns the ordinal number of some
-- basic variable xB[i], which should leave the basis instead the
-- variable y in order to keep primal feasibility. If the adjacent basic
-- solution is primal unbounded and therefore the choice cannot be made,
-- the routine returns zero. */

int lpx_prim_ratio_test(LPX *lp, int len, const int ind[],
      const double val[], int how, double tol)
{     int i, k, m, n, p, t, typx, tagx;
      double alfa_i, abs_alfa_i, big, eps, bbar_i, lb_i, ub_i, temp,
         teta;
      if (!lpx_is_b_avail(lp))
         xfault("lpx_prim_ratio_test: LP basis is not available\n");
      if (lpx_get_prim_stat(lp) != LPX_P_FEAS)
         xfault("lpx_prim_ratio_test: current basic solution is not pri"
            "mal feasible\n");
      if (!(how == +1 || how == -1))
         xfault("lpx_prim_ratio_test: how = %d; invalid parameter\n",
            how);
      m = lpx_get_num_rows(lp);
      n = lpx_get_num_cols(lp);
      /* compute the largest absolute value of the specified influence
         coefficients */
      big = 0.0;
      for (t = 1; t <= len; t++)
      {  temp = val[t];
         if (temp < 0.0) temp = - temp;
         if (big < temp) big = temp;
      }
      /* compute the absolute tolerance eps used to skip small entries
         of the column */
      if (!(0.0 < tol && tol < 1.0))
         xfault("lpx_prim_ratio_test: tol = %g; invalid tolerance\n",
            tol);
      eps = tol * (1.0 + big);
      /* initial settings */
      p = 0, teta = DBL_MAX, big = 0.0;
      /* walk through the entries of the specified column */
      for (t = 1; t <= len; t++)
      {  /* get the ordinal number of basic variable */
         k = ind[t];
         if (!(1 <= k && k <= m+n))
            xfault("lpx_prim_ratio_test: ind[%d] = %d; variable number "
               "out of range\n", t, k);
         if (k <= m)
            tagx = lpx_get_row_stat(lp, k);
         else
            tagx = lpx_get_col_stat(lp, k-m);
         if (tagx != LPX_BS)
            xfault("lpx_prim_ratio_test: ind[%d] = %d; non-basic variab"
               "le not allowed\n", t, k);
         /* determine index of the variable x[k] in the vector xB */
         if (k <= m)
            i = lpx_get_row_b_ind(lp, k);
         else
            i = lpx_get_col_b_ind(lp, k-m);
         xassert(1 <= i && i <= m);
         /* determine unscaled bounds and value of the basic variable
            xB[i] in the current basic solution */
         if (k <= m)
         {  typx = lpx_get_row_type(lp, k);
            lb_i = lpx_get_row_lb(lp, k);
            ub_i = lpx_get_row_ub(lp, k);
            bbar_i = lpx_get_row_prim(lp, k);
         }
         else
         {  typx = lpx_get_col_type(lp, k-m);
            lb_i = lpx_get_col_lb(lp, k-m);
            ub_i = lpx_get_col_ub(lp, k-m);
            bbar_i = lpx_get_col_prim(lp, k-m);
         }
         /* determine influence coefficient for the basic variable
            x[k] = xB[i] in the explicitly specified column and turn to
            the case of increasing the variable y in order to simplify
            the program logic */
         alfa_i = (how > 0 ? +val[t] : -val[t]);
         abs_alfa_i = (alfa_i > 0.0 ? +alfa_i : -alfa_i);
         /* analyze main cases */
         switch (typx)
         {  case LPX_FR:
               /* xB[i] is free variable */
               continue;
            case LPX_LO:
lo:            /* xB[i] has an lower bound */
               if (alfa_i > - eps) continue;
               temp = (lb_i - bbar_i) / alfa_i;
               break;
            case LPX_UP:
up:            /* xB[i] has an upper bound */
               if (alfa_i < + eps) continue;
               temp = (ub_i - bbar_i) / alfa_i;
               break;
            case LPX_DB:
               /* xB[i] has both lower and upper bounds */
               if (alfa_i < 0.0) goto lo; else goto up;
            case LPX_FX:
               /* xB[i] is fixed variable */
               if (abs_alfa_i < eps) continue;
               temp = 0.0;
               break;
            default:
               xassert(typx != typx);
         }
         /* if the value of the variable xB[i] violates its lower or
            upper bound (slightly, because the current basis is assumed
            to be primal feasible), temp is negative; we can think this
            happens due to round-off errors and the value is exactly on
            the bound; this allows replacing temp by zero */
         if (temp < 0.0) temp = 0.0;
         /* apply the minimal ratio test */
         if (teta > temp || teta == temp && big < abs_alfa_i)
            p = k, teta = temp, big = abs_alfa_i;
      }
      /* return the ordinal number of the chosen basic variable */
      return p;
}

/*----------------------------------------------------------------------
-- lpx_dual_ratio_test - perform dual ratio test.
--
-- *Synopsis*
--
-- #include "glplpx.h"
-- int lpx_dual_ratio_test(LPX *lp, int len, int ind[], double val[],
--    int how, double tol);
--
-- *Description*
--
-- The routine lpx_dual_ratio_test performs the dual ratio test for an
-- explicitly specified row of the simplex table.
--
-- The dual basic solution associated with an LP problem object, which
-- the parameter lp points to, should be feasible. No components of the
-- LP problem object are changed by the routine.
--
-- The explicitly specified row of the simplex table is a linear form,
-- which shows how some basic variable y (not necessarily presented in
-- the problem object) depends on non-basic variables xN:
--
--    y = alfa[1]*xN[1] + alfa[2]*xN[2] + ... + alfa[n]*xN[n].       (*)
--
-- The linear form (*) is specified on entry to the routine using the
-- sparse format. Ordinal numbers of non-basic variables xN[j] should be
-- placed in locations ind[1], ..., ind[len], where ordinal numbers 1 to
-- m denote auxiliary variables, and ordinal numbers m+1 to m+n denote
-- structural variables. The corresponding non-zero coefficients alfa[j]
-- should be placed in locations val[1], ..., val[len]. The arrays ind
-- and val are not changed on exit.
--
-- The parameter how specifies in which direction the variable y changes
-- on leaving the basis: +1 means increasing, -1 means decreasing.
--
-- The parameter tol is a relative tolerance (small positive number)
-- used by the routine to skip small alfa[j] of the form (*).
--
-- The routine determines the ordinal number of some non-basic variable
-- (specified in ind[1], ..., ind[len]), which should enter the basis
-- instead the variable y in order to keep dual feasibility, and returns
-- it on exit. If the choice cannot be made (i.e. if the adjacent basic
-- solution is dual unbounded), the routine returns zero.
--
-- *Note*
--
-- If the basic variable y is presented in the LP problem object, the
-- row (*) can be computed using the routine lpx_eval_tab_row. Otherwise
-- it can be computed using the routine lpx_transform_row.
--
-- *Returns*
--
-- The routine lpx_dual_ratio_test returns the ordinal number of some
-- non-basic variable xN[j], which should enter the basis instead the
-- variable y in order to keep dual feasibility. If the adjacent basic
-- solution is dual unbounded and therefore the choice cannot be made,
-- the routine returns zero. */

int lpx_dual_ratio_test(LPX *lp, int len, const int ind[],
      const double val[], int how, double tol)
{     int k, m, n, t, q, tagx;
      double dir, alfa_j, abs_alfa_j, big, eps, cbar_j, temp, teta;
      if (!lpx_is_b_avail(lp))
         xfault("lpx_dual_ratio_test: LP basis is not available\n");
      if (lpx_get_dual_stat(lp) != LPX_D_FEAS)
         xfault("lpx_dual_ratio_test: current basic solution is not dua"
            "l feasible\n");
      if (!(how == +1 || how == -1))
         xfault("lpx_dual_ratio_test: how = %d; invalid parameter\n",
            how);
      m = lpx_get_num_rows(lp);
      n = lpx_get_num_cols(lp);
      dir = (lpx_get_obj_dir(lp) == LPX_MIN ? +1.0 : -1.0);
      /* compute the largest absolute value of the specified influence
         coefficients */
      big = 0.0;
      for (t = 1; t <= len; t++)
      {  temp = val[t];
         if (temp < 0.0) temp = - temp;
         if (big < temp) big = temp;
      }
      /* compute the absolute tolerance eps used to skip small entries
         of the row */
      if (!(0.0 < tol && tol < 1.0))
         xfault("lpx_dual_ratio_test: tol = %g; invalid tolerance\n",
            tol);
      eps = tol * (1.0 + big);
      /* initial settings */
      q = 0, teta = DBL_MAX, big = 0.0;
      /* walk through the entries of the specified row */
      for (t = 1; t <= len; t++)
      {  /* get ordinal number of non-basic variable */
         k = ind[t];
         if (!(1 <= k && k <= m+n))
            xfault("lpx_dual_ratio_test: ind[%d] = %d; variable number "
               "out of range\n", t, k);
         if (k <= m)
            tagx = lpx_get_row_stat(lp, k);
         else
            tagx = lpx_get_col_stat(lp, k-m);
         if (tagx == LPX_BS)
            xfault("lpx_dual_ratio_test: ind[%d] = %d; basic variable n"
               "ot allowed\n", t, k);
         /* determine unscaled reduced cost of the non-basic variable
            x[k] = xN[j] in the current basic solution */
         if (k <= m)
            cbar_j = lpx_get_row_dual(lp, k);
         else
            cbar_j = lpx_get_col_dual(lp, k-m);
         /* determine influence coefficient at the non-basic variable
            x[k] = xN[j] in the explicitly specified row and turn to
            the case of increasing the variable y in order to simplify
            program logic */
         alfa_j = (how > 0 ? +val[t] : -val[t]);
         abs_alfa_j = (alfa_j > 0.0 ? +alfa_j : -alfa_j);
         /* analyze main cases */
         switch (tagx)
         {  case LPX_NL:
               /* xN[j] is on its lower bound */
               if (alfa_j < +eps) continue;
               temp = (dir * cbar_j) / alfa_j;
               break;
            case LPX_NU:
               /* xN[j] is on its upper bound */
               if (alfa_j > -eps) continue;
               temp = (dir * cbar_j) / alfa_j;
               break;
            case LPX_NF:
               /* xN[j] is non-basic free variable */
               if (abs_alfa_j < eps) continue;
               temp = 0.0;
               break;
            case LPX_NS:
               /* xN[j] is non-basic fixed variable */
               continue;
            default:
               xassert(tagx != tagx);
         }
         /* if the reduced cost of the variable xN[j] violates its zero
            bound (slightly, because the current basis is assumed to be
            dual feasible), temp is negative; we can think this happens
            due to round-off errors and the reduced cost is exact zero;
            this allows replacing temp by zero */
         if (temp < 0.0) temp = 0.0;
         /* apply the minimal ratio test */
         if (teta > temp || teta == temp && big < abs_alfa_j)
            q = k, teta = temp, big = abs_alfa_j;
      }
      /* return the ordinal number of the chosen non-basic variable */
      return q;
}

/* eof */