1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/deps/glpk/src/glpnpp01.c	Sun Nov 06 20:59:10 2011 +0100
     1.3 @@ -0,0 +1,927 @@
     1.4 +/* glpnpp01.c */
     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, 2011 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 "glpnpp.h"
    1.29 +
    1.30 +NPP *npp_create_wksp(void)
    1.31 +{     /* create LP/MIP preprocessor workspace */
    1.32 +      NPP *npp;
    1.33 +      npp = xmalloc(sizeof(NPP));
    1.34 +      npp->orig_dir = 0;
    1.35 +      npp->orig_m = npp->orig_n = npp->orig_nnz = 0;
    1.36 +      npp->pool = dmp_create_pool();
    1.37 +      npp->name = npp->obj = NULL;
    1.38 +      npp->c0 = 0.0;
    1.39 +      npp->nrows = npp->ncols = 0;
    1.40 +      npp->r_head = npp->r_tail = NULL;
    1.41 +      npp->c_head = npp->c_tail = NULL;
    1.42 +      npp->stack = dmp_create_pool();
    1.43 +      npp->top = NULL;
    1.44 +#if 0 /* 16/XII-2009 */
    1.45 +      memset(&npp->count, 0, sizeof(npp->count));
    1.46 +#endif
    1.47 +      npp->m = npp->n = npp->nnz = 0;
    1.48 +      npp->row_ref = npp->col_ref = NULL;
    1.49 +      npp->sol = npp->scaling = 0;
    1.50 +      npp->p_stat = npp->d_stat = npp->t_stat = npp->i_stat = 0;
    1.51 +      npp->r_stat = NULL;
    1.52 +      /*npp->r_prim =*/ npp->r_pi = NULL;
    1.53 +      npp->c_stat = NULL;
    1.54 +      npp->c_value = /*npp->c_dual =*/ NULL;
    1.55 +      return npp;
    1.56 +}
    1.57 +
    1.58 +void npp_insert_row(NPP *npp, NPPROW *row, int where)
    1.59 +{     /* insert row to the row list */
    1.60 +      if (where == 0)
    1.61 +      {  /* insert row to the beginning of the row list */
    1.62 +         row->prev = NULL;
    1.63 +         row->next = npp->r_head;
    1.64 +         if (row->next == NULL)
    1.65 +            npp->r_tail = row;
    1.66 +         else
    1.67 +            row->next->prev = row;
    1.68 +         npp->r_head = row;
    1.69 +      }
    1.70 +      else
    1.71 +      {  /* insert row to the end of the row list */
    1.72 +         row->prev = npp->r_tail;
    1.73 +         row->next = NULL;
    1.74 +         if (row->prev == NULL)
    1.75 +            npp->r_head = row;
    1.76 +         else
    1.77 +            row->prev->next = row;
    1.78 +         npp->r_tail = row;
    1.79 +      }
    1.80 +      return;
    1.81 +}
    1.82 +
    1.83 +void npp_remove_row(NPP *npp, NPPROW *row)
    1.84 +{     /* remove row from the row list */
    1.85 +      if (row->prev == NULL)
    1.86 +         npp->r_head = row->next;
    1.87 +      else
    1.88 +         row->prev->next = row->next;
    1.89 +      if (row->next == NULL)
    1.90 +         npp->r_tail = row->prev;
    1.91 +      else
    1.92 +         row->next->prev = row->prev;
    1.93 +      return;
    1.94 +}
    1.95 +
    1.96 +void npp_activate_row(NPP *npp, NPPROW *row)
    1.97 +{     /* make row active */
    1.98 +      if (!row->temp)
    1.99 +      {  row->temp = 1;
   1.100 +         /* move the row to the beginning of the row list */
   1.101 +         npp_remove_row(npp, row);
   1.102 +         npp_insert_row(npp, row, 0);
   1.103 +      }
   1.104 +      return;
   1.105 +}
   1.106 +
   1.107 +void npp_deactivate_row(NPP *npp, NPPROW *row)
   1.108 +{     /* make row inactive */
   1.109 +      if (row->temp)
   1.110 +      {  row->temp = 0;
   1.111 +         /* move the row to the end of the row list */
   1.112 +         npp_remove_row(npp, row);
   1.113 +         npp_insert_row(npp, row, 1);
   1.114 +      }
   1.115 +      return;
   1.116 +}
   1.117 +
   1.118 +void npp_insert_col(NPP *npp, NPPCOL *col, int where)
   1.119 +{     /* insert column to the column list */
   1.120 +      if (where == 0)
   1.121 +      {  /* insert column to the beginning of the column list */
   1.122 +         col->prev = NULL;
   1.123 +         col->next = npp->c_head;
   1.124 +         if (col->next == NULL)
   1.125 +            npp->c_tail = col;
   1.126 +         else
   1.127 +            col->next->prev = col;
   1.128 +         npp->c_head = col;
   1.129 +      }
   1.130 +      else
   1.131 +      {  /* insert column to the end of the column list */
   1.132 +         col->prev = npp->c_tail;
   1.133 +         col->next = NULL;
   1.134 +         if (col->prev == NULL)
   1.135 +            npp->c_head = col;
   1.136 +         else
   1.137 +            col->prev->next = col;
   1.138 +         npp->c_tail = col;
   1.139 +      }
   1.140 +      return;
   1.141 +}
   1.142 +
   1.143 +void npp_remove_col(NPP *npp, NPPCOL *col)
   1.144 +{     /* remove column from the column list */
   1.145 +      if (col->prev == NULL)
   1.146 +         npp->c_head = col->next;
   1.147 +      else
   1.148 +         col->prev->next = col->next;
   1.149 +      if (col->next == NULL)
   1.150 +         npp->c_tail = col->prev;
   1.151 +      else
   1.152 +         col->next->prev = col->prev;
   1.153 +      return;
   1.154 +}
   1.155 +
   1.156 +void npp_activate_col(NPP *npp, NPPCOL *col)
   1.157 +{     /* make column active */
   1.158 +      if (!col->temp)
   1.159 +      {  col->temp = 1;
   1.160 +         /* move the column to the beginning of the column list */
   1.161 +         npp_remove_col(npp, col);
   1.162 +         npp_insert_col(npp, col, 0);
   1.163 +      }
   1.164 +      return;
   1.165 +}
   1.166 +
   1.167 +void npp_deactivate_col(NPP *npp, NPPCOL *col)
   1.168 +{     /* make column inactive */
   1.169 +      if (col->temp)
   1.170 +      {  col->temp = 0;
   1.171 +         /* move the column to the end of the column list */
   1.172 +         npp_remove_col(npp, col);
   1.173 +         npp_insert_col(npp, col, 1);
   1.174 +      }
   1.175 +      return;
   1.176 +}
   1.177 +
   1.178 +NPPROW *npp_add_row(NPP *npp)
   1.179 +{     /* add new row to the current problem */
   1.180 +      NPPROW *row;
   1.181 +      row = dmp_get_atom(npp->pool, sizeof(NPPROW));
   1.182 +      row->i = ++(npp->nrows);
   1.183 +      row->name = NULL;
   1.184 +      row->lb = -DBL_MAX, row->ub = +DBL_MAX;
   1.185 +      row->ptr = NULL;
   1.186 +      row->temp = 0;
   1.187 +      npp_insert_row(npp, row, 1);
   1.188 +      return row;
   1.189 +}
   1.190 +
   1.191 +NPPCOL *npp_add_col(NPP *npp)
   1.192 +{     /* add new column to the current problem */
   1.193 +      NPPCOL *col;
   1.194 +      col = dmp_get_atom(npp->pool, sizeof(NPPCOL));
   1.195 +      col->j = ++(npp->ncols);
   1.196 +      col->name = NULL;
   1.197 +#if 0
   1.198 +      col->kind = GLP_CV;
   1.199 +#else
   1.200 +      col->is_int = 0;
   1.201 +#endif
   1.202 +      col->lb = col->ub = col->coef = 0.0;
   1.203 +      col->ptr = NULL;
   1.204 +      col->temp = 0;
   1.205 +      npp_insert_col(npp, col, 1);
   1.206 +      return col;
   1.207 +}
   1.208 +
   1.209 +NPPAIJ *npp_add_aij(NPP *npp, NPPROW *row, NPPCOL *col, double val)
   1.210 +{     /* add new element to the constraint matrix */
   1.211 +      NPPAIJ *aij;
   1.212 +      aij = dmp_get_atom(npp->pool, sizeof(NPPAIJ));
   1.213 +      aij->row = row;
   1.214 +      aij->col = col;
   1.215 +      aij->val = val;
   1.216 +      aij->r_prev = NULL;
   1.217 +      aij->r_next = row->ptr;
   1.218 +      aij->c_prev = NULL;
   1.219 +      aij->c_next = col->ptr;
   1.220 +      if (aij->r_next != NULL)
   1.221 +         aij->r_next->r_prev = aij;
   1.222 +      if (aij->c_next != NULL)
   1.223 +         aij->c_next->c_prev = aij;
   1.224 +      row->ptr = col->ptr = aij;
   1.225 +      return aij;
   1.226 +}
   1.227 +
   1.228 +int npp_row_nnz(NPP *npp, NPPROW *row)
   1.229 +{     /* count number of non-zero coefficients in row */
   1.230 +      NPPAIJ *aij;
   1.231 +      int nnz;
   1.232 +      xassert(npp == npp);
   1.233 +      nnz = 0;
   1.234 +      for (aij = row->ptr; aij != NULL; aij = aij->r_next)
   1.235 +         nnz++;
   1.236 +      return nnz;
   1.237 +}
   1.238 +
   1.239 +int npp_col_nnz(NPP *npp, NPPCOL *col)
   1.240 +{     /* count number of non-zero coefficients in column */
   1.241 +      NPPAIJ *aij;
   1.242 +      int nnz;
   1.243 +      xassert(npp == npp);
   1.244 +      nnz = 0;
   1.245 +      for (aij = col->ptr; aij != NULL; aij = aij->c_next)
   1.246 +         nnz++;
   1.247 +      return nnz;
   1.248 +}
   1.249 +
   1.250 +void *npp_push_tse(NPP *npp, int (*func)(NPP *npp, void *info),
   1.251 +      int size)
   1.252 +{     /* push new entry to the transformation stack */
   1.253 +      NPPTSE *tse;
   1.254 +      tse = dmp_get_atom(npp->stack, sizeof(NPPTSE));
   1.255 +      tse->func = func;
   1.256 +      tse->info = dmp_get_atom(npp->stack, size);
   1.257 +      tse->link = npp->top;
   1.258 +      npp->top = tse;
   1.259 +      return tse->info;
   1.260 +}
   1.261 +
   1.262 +#if 1 /* 23/XII-2009 */
   1.263 +void npp_erase_row(NPP *npp, NPPROW *row)
   1.264 +{     /* erase row content to make it empty */
   1.265 +      NPPAIJ *aij;
   1.266 +      while (row->ptr != NULL)
   1.267 +      {  aij = row->ptr;
   1.268 +         row->ptr = aij->r_next;
   1.269 +         if (aij->c_prev == NULL)
   1.270 +            aij->col->ptr = aij->c_next;
   1.271 +         else
   1.272 +            aij->c_prev->c_next = aij->c_next;
   1.273 +         if (aij->c_next == NULL)
   1.274 +            ;
   1.275 +         else
   1.276 +            aij->c_next->c_prev = aij->c_prev;
   1.277 +         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
   1.278 +      }
   1.279 +      return;
   1.280 +}
   1.281 +#endif
   1.282 +
   1.283 +void npp_del_row(NPP *npp, NPPROW *row)
   1.284 +{     /* remove row from the current problem */
   1.285 +#if 0 /* 23/XII-2009 */
   1.286 +      NPPAIJ *aij;
   1.287 +#endif
   1.288 +      if (row->name != NULL)
   1.289 +         dmp_free_atom(npp->pool, row->name, strlen(row->name)+1);
   1.290 +#if 0 /* 23/XII-2009 */
   1.291 +      while (row->ptr != NULL)
   1.292 +      {  aij = row->ptr;
   1.293 +         row->ptr = aij->r_next;
   1.294 +         if (aij->c_prev == NULL)
   1.295 +            aij->col->ptr = aij->c_next;
   1.296 +         else
   1.297 +            aij->c_prev->c_next = aij->c_next;
   1.298 +         if (aij->c_next == NULL)
   1.299 +            ;
   1.300 +         else
   1.301 +            aij->c_next->c_prev = aij->c_prev;
   1.302 +         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
   1.303 +      }
   1.304 +#else
   1.305 +      npp_erase_row(npp, row);
   1.306 +#endif
   1.307 +      npp_remove_row(npp, row);
   1.308 +      dmp_free_atom(npp->pool, row, sizeof(NPPROW));
   1.309 +      return;
   1.310 +}
   1.311 +
   1.312 +void npp_del_col(NPP *npp, NPPCOL *col)
   1.313 +{     /* remove column from the current problem */
   1.314 +      NPPAIJ *aij;
   1.315 +      if (col->name != NULL)
   1.316 +         dmp_free_atom(npp->pool, col->name, strlen(col->name)+1);
   1.317 +      while (col->ptr != NULL)
   1.318 +      {  aij = col->ptr;
   1.319 +         col->ptr = aij->c_next;
   1.320 +         if (aij->r_prev == NULL)
   1.321 +            aij->row->ptr = aij->r_next;
   1.322 +         else
   1.323 +            aij->r_prev->r_next = aij->r_next;
   1.324 +         if (aij->r_next == NULL)
   1.325 +            ;
   1.326 +         else
   1.327 +            aij->r_next->r_prev = aij->r_prev;
   1.328 +         dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
   1.329 +      }
   1.330 +      npp_remove_col(npp, col);
   1.331 +      dmp_free_atom(npp->pool, col, sizeof(NPPCOL));
   1.332 +      return;
   1.333 +}
   1.334 +
   1.335 +void npp_del_aij(NPP *npp, NPPAIJ *aij)
   1.336 +{     /* remove element from the constraint matrix */
   1.337 +      if (aij->r_prev == NULL)
   1.338 +         aij->row->ptr = aij->r_next;
   1.339 +      else
   1.340 +         aij->r_prev->r_next = aij->r_next;
   1.341 +      if (aij->r_next == NULL)
   1.342 +         ;
   1.343 +      else
   1.344 +         aij->r_next->r_prev = aij->r_prev;
   1.345 +      if (aij->c_prev == NULL)
   1.346 +         aij->col->ptr = aij->c_next;
   1.347 +      else
   1.348 +         aij->c_prev->c_next = aij->c_next;
   1.349 +      if (aij->c_next == NULL)
   1.350 +         ;
   1.351 +      else
   1.352 +         aij->c_next->c_prev = aij->c_prev;
   1.353 +      dmp_free_atom(npp->pool, aij, sizeof(NPPAIJ));
   1.354 +      return;
   1.355 +}
   1.356 +
   1.357 +void npp_load_prob(NPP *npp, glp_prob *orig, int names, int sol,
   1.358 +      int scaling)
   1.359 +{     /* load original problem into the preprocessor workspace */
   1.360 +      int m = orig->m;
   1.361 +      int n = orig->n;
   1.362 +      NPPROW **link;
   1.363 +      int i, j;
   1.364 +      double dir;
   1.365 +      xassert(names == GLP_OFF || names == GLP_ON);
   1.366 +      xassert(sol == GLP_SOL || sol == GLP_IPT || sol == GLP_MIP);
   1.367 +      xassert(scaling == GLP_OFF || scaling == GLP_ON);
   1.368 +      if (sol == GLP_MIP) xassert(!scaling);
   1.369 +      npp->orig_dir = orig->dir;
   1.370 +      if (npp->orig_dir == GLP_MIN)
   1.371 +         dir = +1.0;
   1.372 +      else if (npp->orig_dir == GLP_MAX)
   1.373 +         dir = -1.0;
   1.374 +      else
   1.375 +         xassert(npp != npp);
   1.376 +      npp->orig_m = m;
   1.377 +      npp->orig_n = n;
   1.378 +      npp->orig_nnz = orig->nnz;
   1.379 +      if (names && orig->name != NULL)
   1.380 +      {  npp->name = dmp_get_atom(npp->pool, strlen(orig->name)+1);
   1.381 +         strcpy(npp->name, orig->name);
   1.382 +      }
   1.383 +      if (names && orig->obj != NULL)
   1.384 +      {  npp->obj = dmp_get_atom(npp->pool, strlen(orig->obj)+1);
   1.385 +         strcpy(npp->obj, orig->obj);
   1.386 +      }
   1.387 +      npp->c0 = dir * orig->c0;
   1.388 +      /* load rows */
   1.389 +      link = xcalloc(1+m, sizeof(NPPROW *));
   1.390 +      for (i = 1; i <= m; i++)
   1.391 +      {  GLPROW *rrr = orig->row[i];
   1.392 +         NPPROW *row;
   1.393 +         link[i] = row = npp_add_row(npp);
   1.394 +         xassert(row->i == i);
   1.395 +         if (names && rrr->name != NULL)
   1.396 +         {  row->name = dmp_get_atom(npp->pool, strlen(rrr->name)+1);
   1.397 +            strcpy(row->name, rrr->name);
   1.398 +         }
   1.399 +         if (!scaling)
   1.400 +         {  if (rrr->type == GLP_FR)
   1.401 +               row->lb = -DBL_MAX, row->ub = +DBL_MAX;
   1.402 +            else if (rrr->type == GLP_LO)
   1.403 +               row->lb = rrr->lb, row->ub = +DBL_MAX;
   1.404 +            else if (rrr->type == GLP_UP)
   1.405 +               row->lb = -DBL_MAX, row->ub = rrr->ub;
   1.406 +            else if (rrr->type == GLP_DB)
   1.407 +               row->lb = rrr->lb, row->ub = rrr->ub;
   1.408 +            else if (rrr->type == GLP_FX)
   1.409 +               row->lb = row->ub = rrr->lb;
   1.410 +            else
   1.411 +               xassert(rrr != rrr);
   1.412 +         }
   1.413 +         else
   1.414 +         {  double rii = rrr->rii;
   1.415 +            if (rrr->type == GLP_FR)
   1.416 +               row->lb = -DBL_MAX, row->ub = +DBL_MAX;
   1.417 +            else if (rrr->type == GLP_LO)
   1.418 +               row->lb = rrr->lb * rii, row->ub = +DBL_MAX;
   1.419 +            else if (rrr->type == GLP_UP)
   1.420 +               row->lb = -DBL_MAX, row->ub = rrr->ub * rii;
   1.421 +            else if (rrr->type == GLP_DB)
   1.422 +               row->lb = rrr->lb * rii, row->ub = rrr->ub * rii;
   1.423 +            else if (rrr->type == GLP_FX)
   1.424 +               row->lb = row->ub = rrr->lb * rii;
   1.425 +            else
   1.426 +               xassert(rrr != rrr);
   1.427 +         }
   1.428 +      }
   1.429 +      /* load columns and constraint coefficients */
   1.430 +      for (j = 1; j <= n; j++)
   1.431 +      {  GLPCOL *ccc = orig->col[j];
   1.432 +         GLPAIJ *aaa;
   1.433 +         NPPCOL *col;
   1.434 +         col = npp_add_col(npp);
   1.435 +         xassert(col->j == j);
   1.436 +         if (names && ccc->name != NULL)
   1.437 +         {  col->name = dmp_get_atom(npp->pool, strlen(ccc->name)+1);
   1.438 +            strcpy(col->name, ccc->name);
   1.439 +         }
   1.440 +         if (sol == GLP_MIP)
   1.441 +#if 0
   1.442 +            col->kind = ccc->kind;
   1.443 +#else
   1.444 +            col->is_int = (char)(ccc->kind == GLP_IV);
   1.445 +#endif
   1.446 +         if (!scaling)
   1.447 +         {  if (ccc->type == GLP_FR)
   1.448 +               col->lb = -DBL_MAX, col->ub = +DBL_MAX;
   1.449 +            else if (ccc->type == GLP_LO)
   1.450 +               col->lb = ccc->lb, col->ub = +DBL_MAX;
   1.451 +            else if (ccc->type == GLP_UP)
   1.452 +               col->lb = -DBL_MAX, col->ub = ccc->ub;
   1.453 +            else if (ccc->type == GLP_DB)
   1.454 +               col->lb = ccc->lb, col->ub = ccc->ub;
   1.455 +            else if (ccc->type == GLP_FX)
   1.456 +               col->lb = col->ub = ccc->lb;
   1.457 +            else
   1.458 +               xassert(ccc != ccc);
   1.459 +            col->coef = dir * ccc->coef;
   1.460 +            for (aaa = ccc->ptr; aaa != NULL; aaa = aaa->c_next)
   1.461 +               npp_add_aij(npp, link[aaa->row->i], col, aaa->val);
   1.462 +         }
   1.463 +         else
   1.464 +         {  double sjj = ccc->sjj;
   1.465 +            if (ccc->type == GLP_FR)
   1.466 +               col->lb = -DBL_MAX, col->ub = +DBL_MAX;
   1.467 +            else if (ccc->type == GLP_LO)
   1.468 +               col->lb = ccc->lb / sjj, col->ub = +DBL_MAX;
   1.469 +            else if (ccc->type == GLP_UP)
   1.470 +               col->lb = -DBL_MAX, col->ub = ccc->ub / sjj;
   1.471 +            else if (ccc->type == GLP_DB)
   1.472 +               col->lb = ccc->lb / sjj, col->ub = ccc->ub / sjj;
   1.473 +            else if (ccc->type == GLP_FX)
   1.474 +               col->lb = col->ub = ccc->lb / sjj;
   1.475 +            else
   1.476 +               xassert(ccc != ccc);
   1.477 +            col->coef = dir * ccc->coef * sjj;
   1.478 +            for (aaa = ccc->ptr; aaa != NULL; aaa = aaa->c_next)
   1.479 +               npp_add_aij(npp, link[aaa->row->i], col,
   1.480 +                  aaa->row->rii * aaa->val * sjj);
   1.481 +         }
   1.482 +      }
   1.483 +      xfree(link);
   1.484 +      /* keep solution indicator and scaling option */
   1.485 +      npp->sol = sol;
   1.486 +      npp->scaling = scaling;
   1.487 +      return;
   1.488 +}
   1.489 +
   1.490 +void npp_build_prob(NPP *npp, glp_prob *prob)
   1.491 +{     /* build resultant (preprocessed) problem */
   1.492 +      NPPROW *row;
   1.493 +      NPPCOL *col;
   1.494 +      NPPAIJ *aij;
   1.495 +      int i, j, type, len, *ind;
   1.496 +      double dir, *val;
   1.497 +      glp_erase_prob(prob);
   1.498 +      glp_set_prob_name(prob, npp->name);
   1.499 +      glp_set_obj_name(prob, npp->obj);
   1.500 +      glp_set_obj_dir(prob, npp->orig_dir);
   1.501 +      if (npp->orig_dir == GLP_MIN)
   1.502 +         dir = +1.0;
   1.503 +      else if (npp->orig_dir == GLP_MAX)
   1.504 +         dir = -1.0;
   1.505 +      else
   1.506 +         xassert(npp != npp);
   1.507 +      glp_set_obj_coef(prob, 0, dir * npp->c0);
   1.508 +      /* build rows */
   1.509 +      for (row = npp->r_head; row != NULL; row = row->next)
   1.510 +      {  row->temp = i = glp_add_rows(prob, 1);
   1.511 +         glp_set_row_name(prob, i, row->name);
   1.512 +         if (row->lb == -DBL_MAX && row->ub == +DBL_MAX)
   1.513 +            type = GLP_FR;
   1.514 +         else if (row->ub == +DBL_MAX)
   1.515 +            type = GLP_LO;
   1.516 +         else if (row->lb == -DBL_MAX)
   1.517 +            type = GLP_UP;
   1.518 +         else if (row->lb != row->ub)
   1.519 +            type = GLP_DB;
   1.520 +         else
   1.521 +            type = GLP_FX;
   1.522 +         glp_set_row_bnds(prob, i, type, row->lb, row->ub);
   1.523 +      }
   1.524 +      /* build columns and the constraint matrix */
   1.525 +      ind = xcalloc(1+prob->m, sizeof(int));
   1.526 +      val = xcalloc(1+prob->m, sizeof(double));
   1.527 +      for (col = npp->c_head; col != NULL; col = col->next)
   1.528 +      {  j = glp_add_cols(prob, 1);
   1.529 +         glp_set_col_name(prob, j, col->name);
   1.530 +#if 0
   1.531 +         glp_set_col_kind(prob, j, col->kind);
   1.532 +#else
   1.533 +         glp_set_col_kind(prob, j, col->is_int ? GLP_IV : GLP_CV);
   1.534 +#endif
   1.535 +         if (col->lb == -DBL_MAX && col->ub == +DBL_MAX)
   1.536 +            type = GLP_FR;
   1.537 +         else if (col->ub == +DBL_MAX)
   1.538 +            type = GLP_LO;
   1.539 +         else if (col->lb == -DBL_MAX)
   1.540 +            type = GLP_UP;
   1.541 +         else if (col->lb != col->ub)
   1.542 +            type = GLP_DB;
   1.543 +         else
   1.544 +            type = GLP_FX;
   1.545 +         glp_set_col_bnds(prob, j, type, col->lb, col->ub);
   1.546 +         glp_set_obj_coef(prob, j, dir * col->coef);
   1.547 +         len = 0;
   1.548 +         for (aij = col->ptr; aij != NULL; aij = aij->c_next)
   1.549 +         {  len++;
   1.550 +            ind[len] = aij->row->temp;
   1.551 +            val[len] = aij->val;
   1.552 +         }
   1.553 +         glp_set_mat_col(prob, j, len, ind, val);
   1.554 +      }
   1.555 +      xfree(ind);
   1.556 +      xfree(val);
   1.557 +      /* resultant problem has been built */
   1.558 +      npp->m = prob->m;
   1.559 +      npp->n = prob->n;
   1.560 +      npp->nnz = prob->nnz;
   1.561 +      npp->row_ref = xcalloc(1+npp->m, sizeof(int));
   1.562 +      npp->col_ref = xcalloc(1+npp->n, sizeof(int));
   1.563 +      for (row = npp->r_head, i = 0; row != NULL; row = row->next)
   1.564 +         npp->row_ref[++i] = row->i;
   1.565 +      for (col = npp->c_head, j = 0; col != NULL; col = col->next)
   1.566 +         npp->col_ref[++j] = col->j;
   1.567 +      /* transformed problem segment is no longer needed */
   1.568 +      dmp_delete_pool(npp->pool), npp->pool = NULL;
   1.569 +      npp->name = npp->obj = NULL;
   1.570 +      npp->c0 = 0.0;
   1.571 +      npp->r_head = npp->r_tail = NULL;
   1.572 +      npp->c_head = npp->c_tail = NULL;
   1.573 +      return;
   1.574 +}
   1.575 +
   1.576 +void npp_postprocess(NPP *npp, glp_prob *prob)
   1.577 +{     /* postprocess solution from the resultant problem */
   1.578 +      GLPROW *row;
   1.579 +      GLPCOL *col;
   1.580 +      NPPTSE *tse;
   1.581 +      int i, j, k;
   1.582 +      double dir;
   1.583 +      xassert(npp->orig_dir == prob->dir);
   1.584 +      if (npp->orig_dir == GLP_MIN)
   1.585 +         dir = +1.0;
   1.586 +      else if (npp->orig_dir == GLP_MAX)
   1.587 +         dir = -1.0;
   1.588 +      else
   1.589 +         xassert(npp != npp);
   1.590 +      xassert(npp->m == prob->m);
   1.591 +      xassert(npp->n == prob->n);
   1.592 +      xassert(npp->nnz == prob->nnz);
   1.593 +      /* copy solution status */
   1.594 +      if (npp->sol == GLP_SOL)
   1.595 +      {  npp->p_stat = prob->pbs_stat;
   1.596 +         npp->d_stat = prob->dbs_stat;
   1.597 +      }
   1.598 +      else if (npp->sol == GLP_IPT)
   1.599 +         npp->t_stat = prob->ipt_stat;
   1.600 +      else if (npp->sol == GLP_MIP)
   1.601 +         npp->i_stat = prob->mip_stat;
   1.602 +      else
   1.603 +         xassert(npp != npp);
   1.604 +      /* allocate solution arrays */
   1.605 +      if (npp->sol == GLP_SOL)
   1.606 +      {  if (npp->r_stat == NULL)
   1.607 +            npp->r_stat = xcalloc(1+npp->nrows, sizeof(char));
   1.608 +         for (i = 1; i <= npp->nrows; i++)
   1.609 +            npp->r_stat[i] = 0;
   1.610 +         if (npp->c_stat == NULL)
   1.611 +            npp->c_stat = xcalloc(1+npp->ncols, sizeof(char));
   1.612 +         for (j = 1; j <= npp->ncols; j++)
   1.613 +            npp->c_stat[j] = 0;
   1.614 +      }
   1.615 +#if 0
   1.616 +      if (npp->r_prim == NULL)
   1.617 +         npp->r_prim = xcalloc(1+npp->nrows, sizeof(double));
   1.618 +      for (i = 1; i <= npp->nrows; i++)
   1.619 +         npp->r_prim[i] = DBL_MAX;
   1.620 +#endif
   1.621 +      if (npp->c_value == NULL)
   1.622 +         npp->c_value = xcalloc(1+npp->ncols, sizeof(double));
   1.623 +      for (j = 1; j <= npp->ncols; j++)
   1.624 +         npp->c_value[j] = DBL_MAX;
   1.625 +      if (npp->sol != GLP_MIP)
   1.626 +      {  if (npp->r_pi == NULL)
   1.627 +            npp->r_pi = xcalloc(1+npp->nrows, sizeof(double));
   1.628 +         for (i = 1; i <= npp->nrows; i++)
   1.629 +            npp->r_pi[i] = DBL_MAX;
   1.630 +#if 0
   1.631 +         if (npp->c_dual == NULL)
   1.632 +            npp->c_dual = xcalloc(1+npp->ncols, sizeof(double));
   1.633 +         for (j = 1; j <= npp->ncols; j++)
   1.634 +            npp->c_dual[j] = DBL_MAX;
   1.635 +#endif
   1.636 +      }
   1.637 +      /* copy solution components from the resultant problem */
   1.638 +      if (npp->sol == GLP_SOL)
   1.639 +      {  for (i = 1; i <= npp->m; i++)
   1.640 +         {  row = prob->row[i];
   1.641 +            k = npp->row_ref[i];
   1.642 +            npp->r_stat[k] = (char)row->stat;
   1.643 +            /*npp->r_prim[k] = row->prim;*/
   1.644 +            npp->r_pi[k] = dir * row->dual;
   1.645 +         }
   1.646 +         for (j = 1; j <= npp->n; j++)
   1.647 +         {  col = prob->col[j];
   1.648 +            k = npp->col_ref[j];
   1.649 +            npp->c_stat[k] = (char)col->stat;
   1.650 +            npp->c_value[k] = col->prim;
   1.651 +            /*npp->c_dual[k] = dir * col->dual;*/
   1.652 +         }
   1.653 +      }
   1.654 +      else if (npp->sol == GLP_IPT)
   1.655 +      {  for (i = 1; i <= npp->m; i++)
   1.656 +         {  row = prob->row[i];
   1.657 +            k = npp->row_ref[i];
   1.658 +            /*npp->r_prim[k] = row->pval;*/
   1.659 +            npp->r_pi[k] = dir * row->dval;
   1.660 +         }
   1.661 +         for (j = 1; j <= npp->n; j++)
   1.662 +         {  col = prob->col[j];
   1.663 +            k = npp->col_ref[j];
   1.664 +            npp->c_value[k] = col->pval;
   1.665 +            /*npp->c_dual[k] = dir * col->dval;*/
   1.666 +         }
   1.667 +      }
   1.668 +      else if (npp->sol == GLP_MIP)
   1.669 +      {
   1.670 +#if 0
   1.671 +         for (i = 1; i <= npp->m; i++)
   1.672 +         {  row = prob->row[i];
   1.673 +            k = npp->row_ref[i];
   1.674 +            /*npp->r_prim[k] = row->mipx;*/
   1.675 +         }
   1.676 +#endif
   1.677 +         for (j = 1; j <= npp->n; j++)
   1.678 +         {  col = prob->col[j];
   1.679 +            k = npp->col_ref[j];
   1.680 +            npp->c_value[k] = col->mipx;
   1.681 +         }
   1.682 +      }
   1.683 +      else
   1.684 +         xassert(npp != npp);
   1.685 +      /* perform postprocessing to construct solution to the original
   1.686 +         problem */
   1.687 +      for (tse = npp->top; tse != NULL; tse = tse->link)
   1.688 +      {  xassert(tse->func != NULL);
   1.689 +         xassert(tse->func(npp, tse->info) == 0);
   1.690 +      }
   1.691 +      return;
   1.692 +}
   1.693 +
   1.694 +void npp_unload_sol(NPP *npp, glp_prob *orig)
   1.695 +{     /* store solution to the original problem */
   1.696 +      GLPROW *row;
   1.697 +      GLPCOL *col;
   1.698 +      int i, j;
   1.699 +      double dir;
   1.700 +      xassert(npp->orig_dir == orig->dir);
   1.701 +      if (npp->orig_dir == GLP_MIN)
   1.702 +         dir = +1.0;
   1.703 +      else if (npp->orig_dir == GLP_MAX)
   1.704 +         dir = -1.0;
   1.705 +      else
   1.706 +         xassert(npp != npp);
   1.707 +      xassert(npp->orig_m == orig->m);
   1.708 +      xassert(npp->orig_n == orig->n);
   1.709 +      xassert(npp->orig_nnz == orig->nnz);
   1.710 +      if (npp->sol == GLP_SOL)
   1.711 +      {  /* store basic solution */
   1.712 +         orig->valid = 0;
   1.713 +         orig->pbs_stat = npp->p_stat;
   1.714 +         orig->dbs_stat = npp->d_stat;
   1.715 +         orig->obj_val = orig->c0;
   1.716 +         orig->some = 0;
   1.717 +         for (i = 1; i <= orig->m; i++)
   1.718 +         {  row = orig->row[i];
   1.719 +            row->stat = npp->r_stat[i];
   1.720 +            if (!npp->scaling)
   1.721 +            {  /*row->prim = npp->r_prim[i];*/
   1.722 +               row->dual = dir * npp->r_pi[i];
   1.723 +            }
   1.724 +            else
   1.725 +            {  /*row->prim = npp->r_prim[i] / row->rii;*/
   1.726 +               row->dual = dir * npp->r_pi[i] * row->rii;
   1.727 +            }
   1.728 +            if (row->stat == GLP_BS)
   1.729 +               row->dual = 0.0;
   1.730 +            else if (row->stat == GLP_NL)
   1.731 +            {  xassert(row->type == GLP_LO || row->type == GLP_DB);
   1.732 +               row->prim = row->lb;
   1.733 +            }
   1.734 +            else if (row->stat == GLP_NU)
   1.735 +            {  xassert(row->type == GLP_UP || row->type == GLP_DB);
   1.736 +               row->prim = row->ub;
   1.737 +            }
   1.738 +            else if (row->stat == GLP_NF)
   1.739 +            {  xassert(row->type == GLP_FR);
   1.740 +               row->prim = 0.0;
   1.741 +            }
   1.742 +            else if (row->stat == GLP_NS)
   1.743 +            {  xassert(row->type == GLP_FX);
   1.744 +               row->prim = row->lb;
   1.745 +            }
   1.746 +            else
   1.747 +               xassert(row != row);
   1.748 +         }
   1.749 +         for (j = 1; j <= orig->n; j++)
   1.750 +         {  col = orig->col[j];
   1.751 +            col->stat = npp->c_stat[j];
   1.752 +            if (!npp->scaling)
   1.753 +            {  col->prim = npp->c_value[j];
   1.754 +               /*col->dual = dir * npp->c_dual[j];*/
   1.755 +            }
   1.756 +            else
   1.757 +            {  col->prim = npp->c_value[j] * col->sjj;
   1.758 +               /*col->dual = dir * npp->c_dual[j] / col->sjj;*/
   1.759 +            }
   1.760 +            if (col->stat == GLP_BS)
   1.761 +               col->dual = 0.0;
   1.762 +#if 1
   1.763 +            else if (col->stat == GLP_NL)
   1.764 +            {  xassert(col->type == GLP_LO || col->type == GLP_DB);
   1.765 +               col->prim = col->lb;
   1.766 +            }
   1.767 +            else if (col->stat == GLP_NU)
   1.768 +            {  xassert(col->type == GLP_UP || col->type == GLP_DB);
   1.769 +               col->prim = col->ub;
   1.770 +            }
   1.771 +            else if (col->stat == GLP_NF)
   1.772 +            {  xassert(col->type == GLP_FR);
   1.773 +               col->prim = 0.0;
   1.774 +            }
   1.775 +            else if (col->stat == GLP_NS)
   1.776 +            {  xassert(col->type == GLP_FX);
   1.777 +               col->prim = col->lb;
   1.778 +            }
   1.779 +            else
   1.780 +               xassert(col != col);
   1.781 +#endif
   1.782 +            orig->obj_val += col->coef * col->prim;
   1.783 +         }
   1.784 +#if 1
   1.785 +         /* compute primal values of inactive rows */
   1.786 +         for (i = 1; i <= orig->m; i++)
   1.787 +         {  row = orig->row[i];
   1.788 +            if (row->stat == GLP_BS)
   1.789 +            {  GLPAIJ *aij;
   1.790 +               double temp;
   1.791 +               temp = 0.0;
   1.792 +               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
   1.793 +                  temp += aij->val * aij->col->prim;
   1.794 +               row->prim = temp;
   1.795 +            }
   1.796 +         }
   1.797 +         /* compute reduced costs of active columns */
   1.798 +         for (j = 1; j <= orig->n; j++)
   1.799 +         {  col = orig->col[j];
   1.800 +            if (col->stat != GLP_BS)
   1.801 +            {  GLPAIJ *aij;
   1.802 +               double temp;
   1.803 +               temp = col->coef;
   1.804 +               for (aij = col->ptr; aij != NULL; aij = aij->c_next)
   1.805 +                  temp -= aij->val * aij->row->dual;
   1.806 +               col->dual = temp;
   1.807 +            }
   1.808 +         }
   1.809 +#endif
   1.810 +      }
   1.811 +      else if (npp->sol == GLP_IPT)
   1.812 +      {  /* store interior-point solution */
   1.813 +         orig->ipt_stat = npp->t_stat;
   1.814 +         orig->ipt_obj = orig->c0;
   1.815 +         for (i = 1; i <= orig->m; i++)
   1.816 +         {  row = orig->row[i];
   1.817 +            if (!npp->scaling)
   1.818 +            {  /*row->pval = npp->r_prim[i];*/
   1.819 +               row->dval = dir * npp->r_pi[i];
   1.820 +            }
   1.821 +            else
   1.822 +            {  /*row->pval = npp->r_prim[i] / row->rii;*/
   1.823 +               row->dval = dir * npp->r_pi[i] * row->rii;
   1.824 +            }
   1.825 +         }
   1.826 +         for (j = 1; j <= orig->n; j++)
   1.827 +         {  col = orig->col[j];
   1.828 +            if (!npp->scaling)
   1.829 +            {  col->pval = npp->c_value[j];
   1.830 +               /*col->dval = dir * npp->c_dual[j];*/
   1.831 +            }
   1.832 +            else
   1.833 +            {  col->pval = npp->c_value[j] * col->sjj;
   1.834 +               /*col->dval = dir * npp->c_dual[j] / col->sjj;*/
   1.835 +            }
   1.836 +            orig->ipt_obj += col->coef * col->pval;
   1.837 +         }
   1.838 +#if 1
   1.839 +         /* compute row primal values */
   1.840 +         for (i = 1; i <= orig->m; i++)
   1.841 +         {  row = orig->row[i];
   1.842 +            {  GLPAIJ *aij;
   1.843 +               double temp;
   1.844 +               temp = 0.0;
   1.845 +               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
   1.846 +                  temp += aij->val * aij->col->pval;
   1.847 +               row->pval = temp;
   1.848 +            }
   1.849 +         }
   1.850 +         /* compute column dual values */
   1.851 +         for (j = 1; j <= orig->n; j++)
   1.852 +         {  col = orig->col[j];
   1.853 +            {  GLPAIJ *aij;
   1.854 +               double temp;
   1.855 +               temp = col->coef;
   1.856 +               for (aij = col->ptr; aij != NULL; aij = aij->c_next)
   1.857 +                  temp -= aij->val * aij->row->dval;
   1.858 +               col->dval = temp;
   1.859 +            }
   1.860 +         }
   1.861 +#endif
   1.862 +      }
   1.863 +      else if (npp->sol == GLP_MIP)
   1.864 +      {  /* store MIP solution */
   1.865 +         xassert(!npp->scaling);
   1.866 +         orig->mip_stat = npp->i_stat;
   1.867 +         orig->mip_obj = orig->c0;
   1.868 +#if 0
   1.869 +         for (i = 1; i <= orig->m; i++)
   1.870 +         {  row = orig->row[i];
   1.871 +            /*row->mipx = npp->r_prim[i];*/
   1.872 +         }
   1.873 +#endif
   1.874 +         for (j = 1; j <= orig->n; j++)
   1.875 +         {  col = orig->col[j];
   1.876 +            col->mipx = npp->c_value[j];
   1.877 +            if (col->kind == GLP_IV)
   1.878 +               xassert(col->mipx == floor(col->mipx));
   1.879 +            orig->mip_obj += col->coef * col->mipx;
   1.880 +         }
   1.881 +#if 1
   1.882 +         /* compute row primal values */
   1.883 +         for (i = 1; i <= orig->m; i++)
   1.884 +         {  row = orig->row[i];
   1.885 +            {  GLPAIJ *aij;
   1.886 +               double temp;
   1.887 +               temp = 0.0;
   1.888 +               for (aij = row->ptr; aij != NULL; aij = aij->r_next)
   1.889 +                  temp += aij->val * aij->col->mipx;
   1.890 +               row->mipx = temp;
   1.891 +            }
   1.892 +         }
   1.893 +#endif
   1.894 +      }
   1.895 +      else
   1.896 +         xassert(npp != npp);
   1.897 +      return;
   1.898 +}
   1.899 +
   1.900 +void npp_delete_wksp(NPP *npp)
   1.901 +{     /* delete LP/MIP preprocessor workspace */
   1.902 +      if (npp->pool != NULL)
   1.903 +         dmp_delete_pool(npp->pool);
   1.904 +      if (npp->stack != NULL)
   1.905 +         dmp_delete_pool(npp->stack);
   1.906 +      if (npp->row_ref != NULL)
   1.907 +         xfree(npp->row_ref);
   1.908 +      if (npp->col_ref != NULL)
   1.909 +         xfree(npp->col_ref);
   1.910 +      if (npp->r_stat != NULL)
   1.911 +         xfree(npp->r_stat);
   1.912 +#if 0
   1.913 +      if (npp->r_prim != NULL)
   1.914 +         xfree(npp->r_prim);
   1.915 +#endif
   1.916 +      if (npp->r_pi != NULL)
   1.917 +         xfree(npp->r_pi);
   1.918 +      if (npp->c_stat != NULL)
   1.919 +         xfree(npp->c_stat);
   1.920 +      if (npp->c_value != NULL)
   1.921 +         xfree(npp->c_value);
   1.922 +#if 0
   1.923 +      if (npp->c_dual != NULL)
   1.924 +         xfree(npp->c_dual);
   1.925 +#endif
   1.926 +      xfree(npp);
   1.927 +      return;
   1.928 +}
   1.929 +
   1.930 +/* eof */