lemon-project-template-glpk

annotate deps/glpk/src/glpios01.c @ 11:4fc6ad2fb8a6

Test GLPK in src/main.cc
author Alpar Juttner <alpar@cs.elte.hu>
date Sun, 06 Nov 2011 21:43:29 +0100
parents
children
rev   line source
alpar@9 1 /* glpios01.c */
alpar@9 2
alpar@9 3 /***********************************************************************
alpar@9 4 * This code is part of GLPK (GNU Linear Programming Kit).
alpar@9 5 *
alpar@9 6 * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
alpar@9 7 * 2009, 2010, 2011 Andrew Makhorin, Department for Applied Informatics,
alpar@9 8 * Moscow Aviation Institute, Moscow, Russia. All rights reserved.
alpar@9 9 * E-mail: <mao@gnu.org>.
alpar@9 10 *
alpar@9 11 * GLPK is free software: you can redistribute it and/or modify it
alpar@9 12 * under the terms of the GNU General Public License as published by
alpar@9 13 * the Free Software Foundation, either version 3 of the License, or
alpar@9 14 * (at your option) any later version.
alpar@9 15 *
alpar@9 16 * GLPK is distributed in the hope that it will be useful, but WITHOUT
alpar@9 17 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
alpar@9 18 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
alpar@9 19 * License for more details.
alpar@9 20 *
alpar@9 21 * You should have received a copy of the GNU General Public License
alpar@9 22 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.
alpar@9 23 ***********************************************************************/
alpar@9 24
alpar@9 25 #include "glpios.h"
alpar@9 26
alpar@9 27 /***********************************************************************
alpar@9 28 * NAME
alpar@9 29 *
alpar@9 30 * ios_create_tree - create branch-and-bound tree
alpar@9 31 *
alpar@9 32 * SYNOPSIS
alpar@9 33 *
alpar@9 34 * #include "glpios.h"
alpar@9 35 * glp_tree *ios_create_tree(glp_prob *mip, const glp_iocp *parm);
alpar@9 36 *
alpar@9 37 * DESCRIPTION
alpar@9 38 *
alpar@9 39 * The routine ios_create_tree creates the branch-and-bound tree.
alpar@9 40 *
alpar@9 41 * Being created the tree consists of the only root subproblem whose
alpar@9 42 * reference number is 1. Note that initially the root subproblem is in
alpar@9 43 * frozen state and therefore needs to be revived.
alpar@9 44 *
alpar@9 45 * RETURNS
alpar@9 46 *
alpar@9 47 * The routine returns a pointer to the tree created. */
alpar@9 48
alpar@9 49 static IOSNPD *new_node(glp_tree *tree, IOSNPD *parent);
alpar@9 50
alpar@9 51 glp_tree *ios_create_tree(glp_prob *mip, const glp_iocp *parm)
alpar@9 52 { int m = mip->m;
alpar@9 53 int n = mip->n;
alpar@9 54 glp_tree *tree;
alpar@9 55 int i, j;
alpar@9 56 xassert(mip->tree == NULL);
alpar@9 57 mip->tree = tree = xmalloc(sizeof(glp_tree));
alpar@9 58 tree->pool = dmp_create_pool();
alpar@9 59 tree->n = n;
alpar@9 60 /* save original problem components */
alpar@9 61 tree->orig_m = m;
alpar@9 62 tree->orig_type = xcalloc(1+m+n, sizeof(char));
alpar@9 63 tree->orig_lb = xcalloc(1+m+n, sizeof(double));
alpar@9 64 tree->orig_ub = xcalloc(1+m+n, sizeof(double));
alpar@9 65 tree->orig_stat = xcalloc(1+m+n, sizeof(char));
alpar@9 66 tree->orig_prim = xcalloc(1+m+n, sizeof(double));
alpar@9 67 tree->orig_dual = xcalloc(1+m+n, sizeof(double));
alpar@9 68 for (i = 1; i <= m; i++)
alpar@9 69 { GLPROW *row = mip->row[i];
alpar@9 70 tree->orig_type[i] = (char)row->type;
alpar@9 71 tree->orig_lb[i] = row->lb;
alpar@9 72 tree->orig_ub[i] = row->ub;
alpar@9 73 tree->orig_stat[i] = (char)row->stat;
alpar@9 74 tree->orig_prim[i] = row->prim;
alpar@9 75 tree->orig_dual[i] = row->dual;
alpar@9 76 }
alpar@9 77 for (j = 1; j <= n; j++)
alpar@9 78 { GLPCOL *col = mip->col[j];
alpar@9 79 tree->orig_type[m+j] = (char)col->type;
alpar@9 80 tree->orig_lb[m+j] = col->lb;
alpar@9 81 tree->orig_ub[m+j] = col->ub;
alpar@9 82 tree->orig_stat[m+j] = (char)col->stat;
alpar@9 83 tree->orig_prim[m+j] = col->prim;
alpar@9 84 tree->orig_dual[m+j] = col->dual;
alpar@9 85 }
alpar@9 86 tree->orig_obj = mip->obj_val;
alpar@9 87 /* initialize the branch-and-bound tree */
alpar@9 88 tree->nslots = 0;
alpar@9 89 tree->avail = 0;
alpar@9 90 tree->slot = NULL;
alpar@9 91 tree->head = tree->tail = NULL;
alpar@9 92 tree->a_cnt = tree->n_cnt = tree->t_cnt = 0;
alpar@9 93 /* the root subproblem is not solved yet, so its final components
alpar@9 94 are unknown so far */
alpar@9 95 tree->root_m = 0;
alpar@9 96 tree->root_type = NULL;
alpar@9 97 tree->root_lb = tree->root_ub = NULL;
alpar@9 98 tree->root_stat = NULL;
alpar@9 99 /* the current subproblem does not exist yet */
alpar@9 100 tree->curr = NULL;
alpar@9 101 tree->mip = mip;
alpar@9 102 /*tree->solved = 0;*/
alpar@9 103 tree->non_int = xcalloc(1+n, sizeof(char));
alpar@9 104 memset(&tree->non_int[1], 0, n);
alpar@9 105 /* arrays to save parent subproblem components will be allocated
alpar@9 106 later */
alpar@9 107 tree->pred_m = tree->pred_max = 0;
alpar@9 108 tree->pred_type = NULL;
alpar@9 109 tree->pred_lb = tree->pred_ub = NULL;
alpar@9 110 tree->pred_stat = NULL;
alpar@9 111 /* cut generator */
alpar@9 112 tree->local = ios_create_pool(tree);
alpar@9 113 /*tree->first_attempt = 1;*/
alpar@9 114 /*tree->max_added_cuts = 0;*/
alpar@9 115 /*tree->min_eff = 0.0;*/
alpar@9 116 /*tree->miss = 0;*/
alpar@9 117 /*tree->just_selected = 0;*/
alpar@9 118 tree->mir_gen = NULL;
alpar@9 119 tree->clq_gen = NULL;
alpar@9 120 /*tree->round = 0;*/
alpar@9 121 #if 0
alpar@9 122 /* create the conflict graph */
alpar@9 123 tree->n_ref = xcalloc(1+n, sizeof(int));
alpar@9 124 memset(&tree->n_ref[1], 0, n * sizeof(int));
alpar@9 125 tree->c_ref = xcalloc(1+n, sizeof(int));
alpar@9 126 memset(&tree->c_ref[1], 0, n * sizeof(int));
alpar@9 127 tree->g = scg_create_graph(0);
alpar@9 128 tree->j_ref = xcalloc(1+tree->g->n_max, sizeof(int));
alpar@9 129 #endif
alpar@9 130 /* pseudocost branching */
alpar@9 131 tree->pcost = NULL;
alpar@9 132 tree->iwrk = xcalloc(1+n, sizeof(int));
alpar@9 133 tree->dwrk = xcalloc(1+n, sizeof(double));
alpar@9 134 /* initialize control parameters */
alpar@9 135 tree->parm = parm;
alpar@9 136 tree->tm_beg = xtime();
alpar@9 137 tree->tm_lag = xlset(0);
alpar@9 138 tree->sol_cnt = 0;
alpar@9 139 /* initialize advanced solver interface */
alpar@9 140 tree->reason = 0;
alpar@9 141 tree->reopt = 0;
alpar@9 142 tree->reinv = 0;
alpar@9 143 tree->br_var = 0;
alpar@9 144 tree->br_sel = 0;
alpar@9 145 tree->child = 0;
alpar@9 146 tree->next_p = 0;
alpar@9 147 /*tree->btrack = NULL;*/
alpar@9 148 tree->stop = 0;
alpar@9 149 /* create the root subproblem, which initially is identical to
alpar@9 150 the original MIP */
alpar@9 151 new_node(tree, NULL);
alpar@9 152 return tree;
alpar@9 153 }
alpar@9 154
alpar@9 155 /***********************************************************************
alpar@9 156 * NAME
alpar@9 157 *
alpar@9 158 * ios_revive_node - revive specified subproblem
alpar@9 159 *
alpar@9 160 * SYNOPSIS
alpar@9 161 *
alpar@9 162 * #include "glpios.h"
alpar@9 163 * void ios_revive_node(glp_tree *tree, int p);
alpar@9 164 *
alpar@9 165 * DESCRIPTION
alpar@9 166 *
alpar@9 167 * The routine ios_revive_node revives the specified subproblem, whose
alpar@9 168 * reference number is p, and thereby makes it the current subproblem.
alpar@9 169 * Note that the specified subproblem must be active. Besides, if the
alpar@9 170 * current subproblem already exists, it must be frozen before reviving
alpar@9 171 * another subproblem. */
alpar@9 172
alpar@9 173 void ios_revive_node(glp_tree *tree, int p)
alpar@9 174 { glp_prob *mip = tree->mip;
alpar@9 175 IOSNPD *node, *root;
alpar@9 176 /* obtain pointer to the specified subproblem */
alpar@9 177 xassert(1 <= p && p <= tree->nslots);
alpar@9 178 node = tree->slot[p].node;
alpar@9 179 xassert(node != NULL);
alpar@9 180 /* the specified subproblem must be active */
alpar@9 181 xassert(node->count == 0);
alpar@9 182 /* the current subproblem must not exist */
alpar@9 183 xassert(tree->curr == NULL);
alpar@9 184 /* the specified subproblem becomes current */
alpar@9 185 tree->curr = node;
alpar@9 186 /*tree->solved = 0;*/
alpar@9 187 /* obtain pointer to the root subproblem */
alpar@9 188 root = tree->slot[1].node;
alpar@9 189 xassert(root != NULL);
alpar@9 190 /* at this point problem object components correspond to the root
alpar@9 191 subproblem, so if the root subproblem should be revived, there
alpar@9 192 is nothing more to do */
alpar@9 193 if (node == root) goto done;
alpar@9 194 xassert(mip->m == tree->root_m);
alpar@9 195 /* build path from the root to the current node */
alpar@9 196 node->temp = NULL;
alpar@9 197 for (node = node; node != NULL; node = node->up)
alpar@9 198 { if (node->up == NULL)
alpar@9 199 xassert(node == root);
alpar@9 200 else
alpar@9 201 node->up->temp = node;
alpar@9 202 }
alpar@9 203 /* go down from the root to the current node and make necessary
alpar@9 204 changes to restore components of the current subproblem */
alpar@9 205 for (node = root; node != NULL; node = node->temp)
alpar@9 206 { int m = mip->m;
alpar@9 207 int n = mip->n;
alpar@9 208 /* if the current node is reached, the problem object at this
alpar@9 209 point corresponds to its parent, so save attributes of rows
alpar@9 210 and columns for the parent subproblem */
alpar@9 211 if (node->temp == NULL)
alpar@9 212 { int i, j;
alpar@9 213 tree->pred_m = m;
alpar@9 214 /* allocate/reallocate arrays, if necessary */
alpar@9 215 if (tree->pred_max < m + n)
alpar@9 216 { int new_size = m + n + 100;
alpar@9 217 if (tree->pred_type != NULL) xfree(tree->pred_type);
alpar@9 218 if (tree->pred_lb != NULL) xfree(tree->pred_lb);
alpar@9 219 if (tree->pred_ub != NULL) xfree(tree->pred_ub);
alpar@9 220 if (tree->pred_stat != NULL) xfree(tree->pred_stat);
alpar@9 221 tree->pred_max = new_size;
alpar@9 222 tree->pred_type = xcalloc(1+new_size, sizeof(char));
alpar@9 223 tree->pred_lb = xcalloc(1+new_size, sizeof(double));
alpar@9 224 tree->pred_ub = xcalloc(1+new_size, sizeof(double));
alpar@9 225 tree->pred_stat = xcalloc(1+new_size, sizeof(char));
alpar@9 226 }
alpar@9 227 /* save row attributes */
alpar@9 228 for (i = 1; i <= m; i++)
alpar@9 229 { GLPROW *row = mip->row[i];
alpar@9 230 tree->pred_type[i] = (char)row->type;
alpar@9 231 tree->pred_lb[i] = row->lb;
alpar@9 232 tree->pred_ub[i] = row->ub;
alpar@9 233 tree->pred_stat[i] = (char)row->stat;
alpar@9 234 }
alpar@9 235 /* save column attributes */
alpar@9 236 for (j = 1; j <= n; j++)
alpar@9 237 { GLPCOL *col = mip->col[j];
alpar@9 238 tree->pred_type[mip->m+j] = (char)col->type;
alpar@9 239 tree->pred_lb[mip->m+j] = col->lb;
alpar@9 240 tree->pred_ub[mip->m+j] = col->ub;
alpar@9 241 tree->pred_stat[mip->m+j] = (char)col->stat;
alpar@9 242 }
alpar@9 243 }
alpar@9 244 /* change bounds of rows and columns */
alpar@9 245 { IOSBND *b;
alpar@9 246 for (b = node->b_ptr; b != NULL; b = b->next)
alpar@9 247 { if (b->k <= m)
alpar@9 248 glp_set_row_bnds(mip, b->k, b->type, b->lb, b->ub);
alpar@9 249 else
alpar@9 250 glp_set_col_bnds(mip, b->k-m, b->type, b->lb, b->ub);
alpar@9 251 }
alpar@9 252 }
alpar@9 253 /* change statuses of rows and columns */
alpar@9 254 { IOSTAT *s;
alpar@9 255 for (s = node->s_ptr; s != NULL; s = s->next)
alpar@9 256 { if (s->k <= m)
alpar@9 257 glp_set_row_stat(mip, s->k, s->stat);
alpar@9 258 else
alpar@9 259 glp_set_col_stat(mip, s->k-m, s->stat);
alpar@9 260 }
alpar@9 261 }
alpar@9 262 /* add new rows */
alpar@9 263 if (node->r_ptr != NULL)
alpar@9 264 { IOSROW *r;
alpar@9 265 IOSAIJ *a;
alpar@9 266 int i, len, *ind;
alpar@9 267 double *val;
alpar@9 268 ind = xcalloc(1+n, sizeof(int));
alpar@9 269 val = xcalloc(1+n, sizeof(double));
alpar@9 270 for (r = node->r_ptr; r != NULL; r = r->next)
alpar@9 271 { i = glp_add_rows(mip, 1);
alpar@9 272 glp_set_row_name(mip, i, r->name);
alpar@9 273 #if 1 /* 20/IX-2008 */
alpar@9 274 xassert(mip->row[i]->level == 0);
alpar@9 275 mip->row[i]->level = node->level;
alpar@9 276 mip->row[i]->origin = r->origin;
alpar@9 277 mip->row[i]->klass = r->klass;
alpar@9 278 #endif
alpar@9 279 glp_set_row_bnds(mip, i, r->type, r->lb, r->ub);
alpar@9 280 len = 0;
alpar@9 281 for (a = r->ptr; a != NULL; a = a->next)
alpar@9 282 len++, ind[len] = a->j, val[len] = a->val;
alpar@9 283 glp_set_mat_row(mip, i, len, ind, val);
alpar@9 284 glp_set_rii(mip, i, r->rii);
alpar@9 285 glp_set_row_stat(mip, i, r->stat);
alpar@9 286 }
alpar@9 287 xfree(ind);
alpar@9 288 xfree(val);
alpar@9 289 }
alpar@9 290 #if 0
alpar@9 291 /* add new edges to the conflict graph */
alpar@9 292 /* add new cliques to the conflict graph */
alpar@9 293 /* (not implemented yet) */
alpar@9 294 xassert(node->own_nn == 0);
alpar@9 295 xassert(node->own_nc == 0);
alpar@9 296 xassert(node->e_ptr == NULL);
alpar@9 297 #endif
alpar@9 298 }
alpar@9 299 /* the specified subproblem has been revived */
alpar@9 300 node = tree->curr;
alpar@9 301 /* delete its bound change list */
alpar@9 302 while (node->b_ptr != NULL)
alpar@9 303 { IOSBND *b;
alpar@9 304 b = node->b_ptr;
alpar@9 305 node->b_ptr = b->next;
alpar@9 306 dmp_free_atom(tree->pool, b, sizeof(IOSBND));
alpar@9 307 }
alpar@9 308 /* delete its status change list */
alpar@9 309 while (node->s_ptr != NULL)
alpar@9 310 { IOSTAT *s;
alpar@9 311 s = node->s_ptr;
alpar@9 312 node->s_ptr = s->next;
alpar@9 313 dmp_free_atom(tree->pool, s, sizeof(IOSTAT));
alpar@9 314 }
alpar@9 315 #if 1 /* 20/XI-2009 */
alpar@9 316 /* delete its row addition list (additional rows may appear, for
alpar@9 317 example, due to branching on GUB constraints */
alpar@9 318 while (node->r_ptr != NULL)
alpar@9 319 { IOSROW *r;
alpar@9 320 r = node->r_ptr;
alpar@9 321 node->r_ptr = r->next;
alpar@9 322 xassert(r->name == NULL);
alpar@9 323 while (r->ptr != NULL)
alpar@9 324 { IOSAIJ *a;
alpar@9 325 a = r->ptr;
alpar@9 326 r->ptr = a->next;
alpar@9 327 dmp_free_atom(tree->pool, a, sizeof(IOSAIJ));
alpar@9 328 }
alpar@9 329 dmp_free_atom(tree->pool, r, sizeof(IOSROW));
alpar@9 330 }
alpar@9 331 #endif
alpar@9 332 done: return;
alpar@9 333 }
alpar@9 334
alpar@9 335 /***********************************************************************
alpar@9 336 * NAME
alpar@9 337 *
alpar@9 338 * ios_freeze_node - freeze current subproblem
alpar@9 339 *
alpar@9 340 * SYNOPSIS
alpar@9 341 *
alpar@9 342 * #include "glpios.h"
alpar@9 343 * void ios_freeze_node(glp_tree *tree);
alpar@9 344 *
alpar@9 345 * DESCRIPTION
alpar@9 346 *
alpar@9 347 * The routine ios_freeze_node freezes the current subproblem. */
alpar@9 348
alpar@9 349 void ios_freeze_node(glp_tree *tree)
alpar@9 350 { glp_prob *mip = tree->mip;
alpar@9 351 int m = mip->m;
alpar@9 352 int n = mip->n;
alpar@9 353 IOSNPD *node;
alpar@9 354 /* obtain pointer to the current subproblem */
alpar@9 355 node = tree->curr;
alpar@9 356 xassert(node != NULL);
alpar@9 357 if (node->up == NULL)
alpar@9 358 { /* freeze the root subproblem */
alpar@9 359 int k;
alpar@9 360 xassert(node->p == 1);
alpar@9 361 xassert(tree->root_m == 0);
alpar@9 362 xassert(tree->root_type == NULL);
alpar@9 363 xassert(tree->root_lb == NULL);
alpar@9 364 xassert(tree->root_ub == NULL);
alpar@9 365 xassert(tree->root_stat == NULL);
alpar@9 366 tree->root_m = m;
alpar@9 367 tree->root_type = xcalloc(1+m+n, sizeof(char));
alpar@9 368 tree->root_lb = xcalloc(1+m+n, sizeof(double));
alpar@9 369 tree->root_ub = xcalloc(1+m+n, sizeof(double));
alpar@9 370 tree->root_stat = xcalloc(1+m+n, sizeof(char));
alpar@9 371 for (k = 1; k <= m+n; k++)
alpar@9 372 { if (k <= m)
alpar@9 373 { GLPROW *row = mip->row[k];
alpar@9 374 tree->root_type[k] = (char)row->type;
alpar@9 375 tree->root_lb[k] = row->lb;
alpar@9 376 tree->root_ub[k] = row->ub;
alpar@9 377 tree->root_stat[k] = (char)row->stat;
alpar@9 378 }
alpar@9 379 else
alpar@9 380 { GLPCOL *col = mip->col[k-m];
alpar@9 381 tree->root_type[k] = (char)col->type;
alpar@9 382 tree->root_lb[k] = col->lb;
alpar@9 383 tree->root_ub[k] = col->ub;
alpar@9 384 tree->root_stat[k] = (char)col->stat;
alpar@9 385 }
alpar@9 386 }
alpar@9 387 }
alpar@9 388 else
alpar@9 389 { /* freeze non-root subproblem */
alpar@9 390 int root_m = tree->root_m;
alpar@9 391 int pred_m = tree->pred_m;
alpar@9 392 int i, j, k;
alpar@9 393 xassert(pred_m <= m);
alpar@9 394 /* build change lists for rows and columns which exist in the
alpar@9 395 parent subproblem */
alpar@9 396 xassert(node->b_ptr == NULL);
alpar@9 397 xassert(node->s_ptr == NULL);
alpar@9 398 for (k = 1; k <= pred_m + n; k++)
alpar@9 399 { int pred_type, pred_stat, type, stat;
alpar@9 400 double pred_lb, pred_ub, lb, ub;
alpar@9 401 /* determine attributes in the parent subproblem */
alpar@9 402 pred_type = tree->pred_type[k];
alpar@9 403 pred_lb = tree->pred_lb[k];
alpar@9 404 pred_ub = tree->pred_ub[k];
alpar@9 405 pred_stat = tree->pred_stat[k];
alpar@9 406 /* determine attributes in the current subproblem */
alpar@9 407 if (k <= pred_m)
alpar@9 408 { GLPROW *row = mip->row[k];
alpar@9 409 type = row->type;
alpar@9 410 lb = row->lb;
alpar@9 411 ub = row->ub;
alpar@9 412 stat = row->stat;
alpar@9 413 }
alpar@9 414 else
alpar@9 415 { GLPCOL *col = mip->col[k - pred_m];
alpar@9 416 type = col->type;
alpar@9 417 lb = col->lb;
alpar@9 418 ub = col->ub;
alpar@9 419 stat = col->stat;
alpar@9 420 }
alpar@9 421 /* save type and bounds of a row/column, if changed */
alpar@9 422 if (!(pred_type == type && pred_lb == lb && pred_ub == ub))
alpar@9 423 { IOSBND *b;
alpar@9 424 b = dmp_get_atom(tree->pool, sizeof(IOSBND));
alpar@9 425 b->k = k;
alpar@9 426 b->type = (unsigned char)type;
alpar@9 427 b->lb = lb;
alpar@9 428 b->ub = ub;
alpar@9 429 b->next = node->b_ptr;
alpar@9 430 node->b_ptr = b;
alpar@9 431 }
alpar@9 432 /* save status of a row/column, if changed */
alpar@9 433 if (pred_stat != stat)
alpar@9 434 { IOSTAT *s;
alpar@9 435 s = dmp_get_atom(tree->pool, sizeof(IOSTAT));
alpar@9 436 s->k = k;
alpar@9 437 s->stat = (unsigned char)stat;
alpar@9 438 s->next = node->s_ptr;
alpar@9 439 node->s_ptr = s;
alpar@9 440 }
alpar@9 441 }
alpar@9 442 /* save new rows added to the current subproblem */
alpar@9 443 xassert(node->r_ptr == NULL);
alpar@9 444 if (pred_m < m)
alpar@9 445 { int i, len, *ind;
alpar@9 446 double *val;
alpar@9 447 ind = xcalloc(1+n, sizeof(int));
alpar@9 448 val = xcalloc(1+n, sizeof(double));
alpar@9 449 for (i = m; i > pred_m; i--)
alpar@9 450 { GLPROW *row = mip->row[i];
alpar@9 451 IOSROW *r;
alpar@9 452 const char *name;
alpar@9 453 r = dmp_get_atom(tree->pool, sizeof(IOSROW));
alpar@9 454 name = glp_get_row_name(mip, i);
alpar@9 455 if (name == NULL)
alpar@9 456 r->name = NULL;
alpar@9 457 else
alpar@9 458 { r->name = dmp_get_atom(tree->pool, strlen(name)+1);
alpar@9 459 strcpy(r->name, name);
alpar@9 460 }
alpar@9 461 #if 1 /* 20/IX-2008 */
alpar@9 462 r->origin = row->origin;
alpar@9 463 r->klass = row->klass;
alpar@9 464 #endif
alpar@9 465 r->type = (unsigned char)row->type;
alpar@9 466 r->lb = row->lb;
alpar@9 467 r->ub = row->ub;
alpar@9 468 r->ptr = NULL;
alpar@9 469 len = glp_get_mat_row(mip, i, ind, val);
alpar@9 470 for (k = 1; k <= len; k++)
alpar@9 471 { IOSAIJ *a;
alpar@9 472 a = dmp_get_atom(tree->pool, sizeof(IOSAIJ));
alpar@9 473 a->j = ind[k];
alpar@9 474 a->val = val[k];
alpar@9 475 a->next = r->ptr;
alpar@9 476 r->ptr = a;
alpar@9 477 }
alpar@9 478 r->rii = row->rii;
alpar@9 479 r->stat = (unsigned char)row->stat;
alpar@9 480 r->next = node->r_ptr;
alpar@9 481 node->r_ptr = r;
alpar@9 482 }
alpar@9 483 xfree(ind);
alpar@9 484 xfree(val);
alpar@9 485 }
alpar@9 486 /* remove all rows missing in the root subproblem */
alpar@9 487 if (m != root_m)
alpar@9 488 { int nrs, *num;
alpar@9 489 nrs = m - root_m;
alpar@9 490 xassert(nrs > 0);
alpar@9 491 num = xcalloc(1+nrs, sizeof(int));
alpar@9 492 for (i = 1; i <= nrs; i++) num[i] = root_m + i;
alpar@9 493 glp_del_rows(mip, nrs, num);
alpar@9 494 xfree(num);
alpar@9 495 }
alpar@9 496 m = mip->m;
alpar@9 497 /* and restore attributes of all rows and columns for the root
alpar@9 498 subproblem */
alpar@9 499 xassert(m == root_m);
alpar@9 500 for (i = 1; i <= m; i++)
alpar@9 501 { glp_set_row_bnds(mip, i, tree->root_type[i],
alpar@9 502 tree->root_lb[i], tree->root_ub[i]);
alpar@9 503 glp_set_row_stat(mip, i, tree->root_stat[i]);
alpar@9 504 }
alpar@9 505 for (j = 1; j <= n; j++)
alpar@9 506 { glp_set_col_bnds(mip, j, tree->root_type[m+j],
alpar@9 507 tree->root_lb[m+j], tree->root_ub[m+j]);
alpar@9 508 glp_set_col_stat(mip, j, tree->root_stat[m+j]);
alpar@9 509 }
alpar@9 510 #if 1
alpar@9 511 /* remove all edges and cliques missing in the conflict graph
alpar@9 512 for the root subproblem */
alpar@9 513 /* (not implemented yet) */
alpar@9 514 #endif
alpar@9 515 }
alpar@9 516 /* the current subproblem has been frozen */
alpar@9 517 tree->curr = NULL;
alpar@9 518 return;
alpar@9 519 }
alpar@9 520
alpar@9 521 /***********************************************************************
alpar@9 522 * NAME
alpar@9 523 *
alpar@9 524 * ios_clone_node - clone specified subproblem
alpar@9 525 *
alpar@9 526 * SYNOPSIS
alpar@9 527 *
alpar@9 528 * #include "glpios.h"
alpar@9 529 * void ios_clone_node(glp_tree *tree, int p, int nnn, int ref[]);
alpar@9 530 *
alpar@9 531 * DESCRIPTION
alpar@9 532 *
alpar@9 533 * The routine ios_clone_node clones the specified subproblem, whose
alpar@9 534 * reference number is p, creating its nnn exact copies. Note that the
alpar@9 535 * specified subproblem must be active and must be in the frozen state
alpar@9 536 * (i.e. it must not be the current subproblem).
alpar@9 537 *
alpar@9 538 * Each clone, an exact copy of the specified subproblem, becomes a new
alpar@9 539 * active subproblem added to the end of the active list. After cloning
alpar@9 540 * the specified subproblem becomes inactive.
alpar@9 541 *
alpar@9 542 * The reference numbers of clone subproblems are stored to locations
alpar@9 543 * ref[1], ..., ref[nnn]. */
alpar@9 544
alpar@9 545 static int get_slot(glp_tree *tree)
alpar@9 546 { int p;
alpar@9 547 /* if no free slots are available, increase the room */
alpar@9 548 if (tree->avail == 0)
alpar@9 549 { int nslots = tree->nslots;
alpar@9 550 IOSLOT *save = tree->slot;
alpar@9 551 if (nslots == 0)
alpar@9 552 tree->nslots = 20;
alpar@9 553 else
alpar@9 554 { tree->nslots = nslots + nslots;
alpar@9 555 xassert(tree->nslots > nslots);
alpar@9 556 }
alpar@9 557 tree->slot = xcalloc(1+tree->nslots, sizeof(IOSLOT));
alpar@9 558 if (save != NULL)
alpar@9 559 { memcpy(&tree->slot[1], &save[1], nslots * sizeof(IOSLOT));
alpar@9 560 xfree(save);
alpar@9 561 }
alpar@9 562 /* push more free slots into the stack */
alpar@9 563 for (p = tree->nslots; p > nslots; p--)
alpar@9 564 { tree->slot[p].node = NULL;
alpar@9 565 tree->slot[p].next = tree->avail;
alpar@9 566 tree->avail = p;
alpar@9 567 }
alpar@9 568 }
alpar@9 569 /* pull a free slot from the stack */
alpar@9 570 p = tree->avail;
alpar@9 571 tree->avail = tree->slot[p].next;
alpar@9 572 xassert(tree->slot[p].node == NULL);
alpar@9 573 tree->slot[p].next = 0;
alpar@9 574 return p;
alpar@9 575 }
alpar@9 576
alpar@9 577 static IOSNPD *new_node(glp_tree *tree, IOSNPD *parent)
alpar@9 578 { IOSNPD *node;
alpar@9 579 int p;
alpar@9 580 /* pull a free slot for the new node */
alpar@9 581 p = get_slot(tree);
alpar@9 582 /* create descriptor of the new subproblem */
alpar@9 583 node = dmp_get_atom(tree->pool, sizeof(IOSNPD));
alpar@9 584 tree->slot[p].node = node;
alpar@9 585 node->p = p;
alpar@9 586 node->up = parent;
alpar@9 587 node->level = (parent == NULL ? 0 : parent->level + 1);
alpar@9 588 node->count = 0;
alpar@9 589 node->b_ptr = NULL;
alpar@9 590 node->s_ptr = NULL;
alpar@9 591 node->r_ptr = NULL;
alpar@9 592 node->solved = 0;
alpar@9 593 #if 0
alpar@9 594 node->own_nn = node->own_nc = 0;
alpar@9 595 node->e_ptr = NULL;
alpar@9 596 #endif
alpar@9 597 #if 1 /* 04/X-2008 */
alpar@9 598 node->lp_obj = (parent == NULL ? (tree->mip->dir == GLP_MIN ?
alpar@9 599 -DBL_MAX : +DBL_MAX) : parent->lp_obj);
alpar@9 600 #endif
alpar@9 601 node->bound = (parent == NULL ? (tree->mip->dir == GLP_MIN ?
alpar@9 602 -DBL_MAX : +DBL_MAX) : parent->bound);
alpar@9 603 node->br_var = 0;
alpar@9 604 node->br_val = 0.0;
alpar@9 605 node->ii_cnt = 0;
alpar@9 606 node->ii_sum = 0.0;
alpar@9 607 #if 1 /* 30/XI-2009 */
alpar@9 608 node->changed = 0;
alpar@9 609 #endif
alpar@9 610 if (tree->parm->cb_size == 0)
alpar@9 611 node->data = NULL;
alpar@9 612 else
alpar@9 613 { node->data = dmp_get_atom(tree->pool, tree->parm->cb_size);
alpar@9 614 memset(node->data, 0, tree->parm->cb_size);
alpar@9 615 }
alpar@9 616 node->temp = NULL;
alpar@9 617 node->prev = tree->tail;
alpar@9 618 node->next = NULL;
alpar@9 619 /* add the new subproblem to the end of the active list */
alpar@9 620 if (tree->head == NULL)
alpar@9 621 tree->head = node;
alpar@9 622 else
alpar@9 623 tree->tail->next = node;
alpar@9 624 tree->tail = node;
alpar@9 625 tree->a_cnt++;
alpar@9 626 tree->n_cnt++;
alpar@9 627 tree->t_cnt++;
alpar@9 628 /* increase the number of child subproblems */
alpar@9 629 if (parent == NULL)
alpar@9 630 xassert(p == 1);
alpar@9 631 else
alpar@9 632 parent->count++;
alpar@9 633 return node;
alpar@9 634 }
alpar@9 635
alpar@9 636 void ios_clone_node(glp_tree *tree, int p, int nnn, int ref[])
alpar@9 637 { IOSNPD *node;
alpar@9 638 int k;
alpar@9 639 /* obtain pointer to the subproblem to be cloned */
alpar@9 640 xassert(1 <= p && p <= tree->nslots);
alpar@9 641 node = tree->slot[p].node;
alpar@9 642 xassert(node != NULL);
alpar@9 643 /* the specified subproblem must be active */
alpar@9 644 xassert(node->count == 0);
alpar@9 645 /* and must be in the frozen state */
alpar@9 646 xassert(tree->curr != node);
alpar@9 647 /* remove the specified subproblem from the active list, because
alpar@9 648 it becomes inactive */
alpar@9 649 if (node->prev == NULL)
alpar@9 650 tree->head = node->next;
alpar@9 651 else
alpar@9 652 node->prev->next = node->next;
alpar@9 653 if (node->next == NULL)
alpar@9 654 tree->tail = node->prev;
alpar@9 655 else
alpar@9 656 node->next->prev = node->prev;
alpar@9 657 node->prev = node->next = NULL;
alpar@9 658 tree->a_cnt--;
alpar@9 659 /* create clone subproblems */
alpar@9 660 xassert(nnn > 0);
alpar@9 661 for (k = 1; k <= nnn; k++)
alpar@9 662 ref[k] = new_node(tree, node)->p;
alpar@9 663 return;
alpar@9 664 }
alpar@9 665
alpar@9 666 /***********************************************************************
alpar@9 667 * NAME
alpar@9 668 *
alpar@9 669 * ios_delete_node - delete specified subproblem
alpar@9 670 *
alpar@9 671 * SYNOPSIS
alpar@9 672 *
alpar@9 673 * #include "glpios.h"
alpar@9 674 * void ios_delete_node(glp_tree *tree, int p);
alpar@9 675 *
alpar@9 676 * DESCRIPTION
alpar@9 677 *
alpar@9 678 * The routine ios_delete_node deletes the specified subproblem, whose
alpar@9 679 * reference number is p. The subproblem must be active and must be in
alpar@9 680 * the frozen state (i.e. it must not be the current subproblem).
alpar@9 681 *
alpar@9 682 * Note that deletion is performed recursively, i.e. if a subproblem to
alpar@9 683 * be deleted is the only child of its parent, the parent subproblem is
alpar@9 684 * also deleted, etc. */
alpar@9 685
alpar@9 686 void ios_delete_node(glp_tree *tree, int p)
alpar@9 687 { IOSNPD *node, *temp;
alpar@9 688 /* obtain pointer to the subproblem to be deleted */
alpar@9 689 xassert(1 <= p && p <= tree->nslots);
alpar@9 690 node = tree->slot[p].node;
alpar@9 691 xassert(node != NULL);
alpar@9 692 /* the specified subproblem must be active */
alpar@9 693 xassert(node->count == 0);
alpar@9 694 /* and must be in the frozen state */
alpar@9 695 xassert(tree->curr != node);
alpar@9 696 /* remove the specified subproblem from the active list, because
alpar@9 697 it is gone from the tree */
alpar@9 698 if (node->prev == NULL)
alpar@9 699 tree->head = node->next;
alpar@9 700 else
alpar@9 701 node->prev->next = node->next;
alpar@9 702 if (node->next == NULL)
alpar@9 703 tree->tail = node->prev;
alpar@9 704 else
alpar@9 705 node->next->prev = node->prev;
alpar@9 706 node->prev = node->next = NULL;
alpar@9 707 tree->a_cnt--;
alpar@9 708 loop: /* recursive deletion starts here */
alpar@9 709 /* delete the bound change list */
alpar@9 710 { IOSBND *b;
alpar@9 711 while (node->b_ptr != NULL)
alpar@9 712 { b = node->b_ptr;
alpar@9 713 node->b_ptr = b->next;
alpar@9 714 dmp_free_atom(tree->pool, b, sizeof(IOSBND));
alpar@9 715 }
alpar@9 716 }
alpar@9 717 /* delete the status change list */
alpar@9 718 { IOSTAT *s;
alpar@9 719 while (node->s_ptr != NULL)
alpar@9 720 { s = node->s_ptr;
alpar@9 721 node->s_ptr = s->next;
alpar@9 722 dmp_free_atom(tree->pool, s, sizeof(IOSTAT));
alpar@9 723 }
alpar@9 724 }
alpar@9 725 /* delete the row addition list */
alpar@9 726 while (node->r_ptr != NULL)
alpar@9 727 { IOSROW *r;
alpar@9 728 r = node->r_ptr;
alpar@9 729 if (r->name != NULL)
alpar@9 730 dmp_free_atom(tree->pool, r->name, strlen(r->name)+1);
alpar@9 731 while (r->ptr != NULL)
alpar@9 732 { IOSAIJ *a;
alpar@9 733 a = r->ptr;
alpar@9 734 r->ptr = a->next;
alpar@9 735 dmp_free_atom(tree->pool, a, sizeof(IOSAIJ));
alpar@9 736 }
alpar@9 737 node->r_ptr = r->next;
alpar@9 738 dmp_free_atom(tree->pool, r, sizeof(IOSROW));
alpar@9 739 }
alpar@9 740 #if 0
alpar@9 741 /* delete the edge addition list */
alpar@9 742 /* delete the clique addition list */
alpar@9 743 /* (not implemented yet) */
alpar@9 744 xassert(node->own_nn == 0);
alpar@9 745 xassert(node->own_nc == 0);
alpar@9 746 xassert(node->e_ptr == NULL);
alpar@9 747 #endif
alpar@9 748 /* free application-specific data */
alpar@9 749 if (tree->parm->cb_size == 0)
alpar@9 750 xassert(node->data == NULL);
alpar@9 751 else
alpar@9 752 dmp_free_atom(tree->pool, node->data, tree->parm->cb_size);
alpar@9 753 /* free the corresponding node slot */
alpar@9 754 p = node->p;
alpar@9 755 xassert(tree->slot[p].node == node);
alpar@9 756 tree->slot[p].node = NULL;
alpar@9 757 tree->slot[p].next = tree->avail;
alpar@9 758 tree->avail = p;
alpar@9 759 /* save pointer to the parent subproblem */
alpar@9 760 temp = node->up;
alpar@9 761 /* delete the subproblem descriptor */
alpar@9 762 dmp_free_atom(tree->pool, node, sizeof(IOSNPD));
alpar@9 763 tree->n_cnt--;
alpar@9 764 /* take pointer to the parent subproblem */
alpar@9 765 node = temp;
alpar@9 766 if (node != NULL)
alpar@9 767 { /* the parent subproblem exists; decrease the number of its
alpar@9 768 child subproblems */
alpar@9 769 xassert(node->count > 0);
alpar@9 770 node->count--;
alpar@9 771 /* if now the parent subproblem has no childs, it also must be
alpar@9 772 deleted */
alpar@9 773 if (node->count == 0) goto loop;
alpar@9 774 }
alpar@9 775 return;
alpar@9 776 }
alpar@9 777
alpar@9 778 /***********************************************************************
alpar@9 779 * NAME
alpar@9 780 *
alpar@9 781 * ios_delete_tree - delete branch-and-bound tree
alpar@9 782 *
alpar@9 783 * SYNOPSIS
alpar@9 784 *
alpar@9 785 * #include "glpios.h"
alpar@9 786 * void ios_delete_tree(glp_tree *tree);
alpar@9 787 *
alpar@9 788 * DESCRIPTION
alpar@9 789 *
alpar@9 790 * The routine ios_delete_tree deletes the branch-and-bound tree, which
alpar@9 791 * the parameter tree points to, and frees all the memory allocated to
alpar@9 792 * this program object.
alpar@9 793 *
alpar@9 794 * On exit components of the problem object are restored to correspond
alpar@9 795 * to the original MIP passed to the routine ios_create_tree. */
alpar@9 796
alpar@9 797 void ios_delete_tree(glp_tree *tree)
alpar@9 798 { glp_prob *mip = tree->mip;
alpar@9 799 int i, j;
alpar@9 800 int m = mip->m;
alpar@9 801 int n = mip->n;
alpar@9 802 xassert(mip->tree == tree);
alpar@9 803 /* remove all additional rows */
alpar@9 804 if (m != tree->orig_m)
alpar@9 805 { int nrs, *num;
alpar@9 806 nrs = m - tree->orig_m;
alpar@9 807 xassert(nrs > 0);
alpar@9 808 num = xcalloc(1+nrs, sizeof(int));
alpar@9 809 for (i = 1; i <= nrs; i++) num[i] = tree->orig_m + i;
alpar@9 810 glp_del_rows(mip, nrs, num);
alpar@9 811 xfree(num);
alpar@9 812 }
alpar@9 813 m = tree->orig_m;
alpar@9 814 /* restore original attributes of rows and columns */
alpar@9 815 xassert(m == tree->orig_m);
alpar@9 816 xassert(n == tree->n);
alpar@9 817 for (i = 1; i <= m; i++)
alpar@9 818 { glp_set_row_bnds(mip, i, tree->orig_type[i],
alpar@9 819 tree->orig_lb[i], tree->orig_ub[i]);
alpar@9 820 glp_set_row_stat(mip, i, tree->orig_stat[i]);
alpar@9 821 mip->row[i]->prim = tree->orig_prim[i];
alpar@9 822 mip->row[i]->dual = tree->orig_dual[i];
alpar@9 823 }
alpar@9 824 for (j = 1; j <= n; j++)
alpar@9 825 { glp_set_col_bnds(mip, j, tree->orig_type[m+j],
alpar@9 826 tree->orig_lb[m+j], tree->orig_ub[m+j]);
alpar@9 827 glp_set_col_stat(mip, j, tree->orig_stat[m+j]);
alpar@9 828 mip->col[j]->prim = tree->orig_prim[m+j];
alpar@9 829 mip->col[j]->dual = tree->orig_dual[m+j];
alpar@9 830 }
alpar@9 831 mip->pbs_stat = mip->dbs_stat = GLP_FEAS;
alpar@9 832 mip->obj_val = tree->orig_obj;
alpar@9 833 /* delete the branch-and-bound tree */
alpar@9 834 xassert(tree->local != NULL);
alpar@9 835 ios_delete_pool(tree, tree->local);
alpar@9 836 dmp_delete_pool(tree->pool);
alpar@9 837 xfree(tree->orig_type);
alpar@9 838 xfree(tree->orig_lb);
alpar@9 839 xfree(tree->orig_ub);
alpar@9 840 xfree(tree->orig_stat);
alpar@9 841 xfree(tree->orig_prim);
alpar@9 842 xfree(tree->orig_dual);
alpar@9 843 xfree(tree->slot);
alpar@9 844 if (tree->root_type != NULL) xfree(tree->root_type);
alpar@9 845 if (tree->root_lb != NULL) xfree(tree->root_lb);
alpar@9 846 if (tree->root_ub != NULL) xfree(tree->root_ub);
alpar@9 847 if (tree->root_stat != NULL) xfree(tree->root_stat);
alpar@9 848 xfree(tree->non_int);
alpar@9 849 #if 0
alpar@9 850 xfree(tree->n_ref);
alpar@9 851 xfree(tree->c_ref);
alpar@9 852 xfree(tree->j_ref);
alpar@9 853 #endif
alpar@9 854 if (tree->pcost != NULL) ios_pcost_free(tree);
alpar@9 855 xfree(tree->iwrk);
alpar@9 856 xfree(tree->dwrk);
alpar@9 857 #if 0
alpar@9 858 scg_delete_graph(tree->g);
alpar@9 859 #endif
alpar@9 860 if (tree->pred_type != NULL) xfree(tree->pred_type);
alpar@9 861 if (tree->pred_lb != NULL) xfree(tree->pred_lb);
alpar@9 862 if (tree->pred_ub != NULL) xfree(tree->pred_ub);
alpar@9 863 if (tree->pred_stat != NULL) xfree(tree->pred_stat);
alpar@9 864 #if 0
alpar@9 865 xassert(tree->cut_gen == NULL);
alpar@9 866 #endif
alpar@9 867 xassert(tree->mir_gen == NULL);
alpar@9 868 xassert(tree->clq_gen == NULL);
alpar@9 869 xfree(tree);
alpar@9 870 mip->tree = NULL;
alpar@9 871 return;
alpar@9 872 }
alpar@9 873
alpar@9 874 /***********************************************************************
alpar@9 875 * NAME
alpar@9 876 *
alpar@9 877 * ios_eval_degrad - estimate obj. degrad. for down- and up-branches
alpar@9 878 *
alpar@9 879 * SYNOPSIS
alpar@9 880 *
alpar@9 881 * #include "glpios.h"
alpar@9 882 * void ios_eval_degrad(glp_tree *tree, int j, double *dn, double *up);
alpar@9 883 *
alpar@9 884 * DESCRIPTION
alpar@9 885 *
alpar@9 886 * Given optimal basis to LP relaxation of the current subproblem the
alpar@9 887 * routine ios_eval_degrad performs the dual ratio test to compute the
alpar@9 888 * objective values in the adjacent basis for down- and up-branches,
alpar@9 889 * which are stored in locations *dn and *up, assuming that x[j] is a
alpar@9 890 * variable chosen to branch upon. */
alpar@9 891
alpar@9 892 void ios_eval_degrad(glp_tree *tree, int j, double *dn, double *up)
alpar@9 893 { glp_prob *mip = tree->mip;
alpar@9 894 int m = mip->m, n = mip->n;
alpar@9 895 int len, kase, k, t, stat;
alpar@9 896 double alfa, beta, gamma, delta, dz;
alpar@9 897 int *ind = tree->iwrk;
alpar@9 898 double *val = tree->dwrk;
alpar@9 899 /* current basis must be optimal */
alpar@9 900 xassert(glp_get_status(mip) == GLP_OPT);
alpar@9 901 /* basis factorization must exist */
alpar@9 902 xassert(glp_bf_exists(mip));
alpar@9 903 /* obtain (fractional) value of x[j] in optimal basic solution
alpar@9 904 to LP relaxation of the current subproblem */
alpar@9 905 xassert(1 <= j && j <= n);
alpar@9 906 beta = mip->col[j]->prim;
alpar@9 907 /* since the value of x[j] is fractional, it is basic; compute
alpar@9 908 corresponding row of the simplex table */
alpar@9 909 len = lpx_eval_tab_row(mip, m+j, ind, val);
alpar@9 910 /* kase < 0 means down-branch; kase > 0 means up-branch */
alpar@9 911 for (kase = -1; kase <= +1; kase += 2)
alpar@9 912 { /* for down-branch we introduce new upper bound floor(beta)
alpar@9 913 for x[j]; similarly, for up-branch we introduce new lower
alpar@9 914 bound ceil(beta) for x[j]; in the current basis this new
alpar@9 915 upper/lower bound is violated, so in the adjacent basis
alpar@9 916 x[j] will leave the basis and go to its new upper/lower
alpar@9 917 bound; we need to know which non-basic variable x[k] should
alpar@9 918 enter the basis to keep dual feasibility */
alpar@9 919 #if 0 /* 23/XI-2009 */
alpar@9 920 k = lpx_dual_ratio_test(mip, len, ind, val, kase, 1e-7);
alpar@9 921 #else
alpar@9 922 k = lpx_dual_ratio_test(mip, len, ind, val, kase, 1e-9);
alpar@9 923 #endif
alpar@9 924 /* if no variable has been chosen, current basis being primal
alpar@9 925 infeasible due to the new upper/lower bound of x[j] is dual
alpar@9 926 unbounded, therefore, LP relaxation to corresponding branch
alpar@9 927 has no primal feasible solution */
alpar@9 928 if (k == 0)
alpar@9 929 { if (mip->dir == GLP_MIN)
alpar@9 930 { if (kase < 0)
alpar@9 931 *dn = +DBL_MAX;
alpar@9 932 else
alpar@9 933 *up = +DBL_MAX;
alpar@9 934 }
alpar@9 935 else if (mip->dir == GLP_MAX)
alpar@9 936 { if (kase < 0)
alpar@9 937 *dn = -DBL_MAX;
alpar@9 938 else
alpar@9 939 *up = -DBL_MAX;
alpar@9 940 }
alpar@9 941 else
alpar@9 942 xassert(mip != mip);
alpar@9 943 continue;
alpar@9 944 }
alpar@9 945 xassert(1 <= k && k <= m+n);
alpar@9 946 /* row of the simplex table corresponding to specified basic
alpar@9 947 variable x[j] is the following:
alpar@9 948 x[j] = ... + alfa * x[k] + ... ;
alpar@9 949 we need to know influence coefficient, alfa, at non-basic
alpar@9 950 variable x[k] chosen with the dual ratio test */
alpar@9 951 for (t = 1; t <= len; t++)
alpar@9 952 if (ind[t] == k) break;
alpar@9 953 xassert(1 <= t && t <= len);
alpar@9 954 alfa = val[t];
alpar@9 955 /* determine status and reduced cost of variable x[k] */
alpar@9 956 if (k <= m)
alpar@9 957 { stat = mip->row[k]->stat;
alpar@9 958 gamma = mip->row[k]->dual;
alpar@9 959 }
alpar@9 960 else
alpar@9 961 { stat = mip->col[k-m]->stat;
alpar@9 962 gamma = mip->col[k-m]->dual;
alpar@9 963 }
alpar@9 964 /* x[k] cannot be basic or fixed non-basic */
alpar@9 965 xassert(stat == GLP_NL || stat == GLP_NU || stat == GLP_NF);
alpar@9 966 /* if the current basis is dual degenerative, some reduced
alpar@9 967 costs, which are close to zero, may have wrong sign due to
alpar@9 968 round-off errors, so correct the sign of gamma */
alpar@9 969 if (mip->dir == GLP_MIN)
alpar@9 970 { if (stat == GLP_NL && gamma < 0.0 ||
alpar@9 971 stat == GLP_NU && gamma > 0.0 ||
alpar@9 972 stat == GLP_NF) gamma = 0.0;
alpar@9 973 }
alpar@9 974 else if (mip->dir == GLP_MAX)
alpar@9 975 { if (stat == GLP_NL && gamma > 0.0 ||
alpar@9 976 stat == GLP_NU && gamma < 0.0 ||
alpar@9 977 stat == GLP_NF) gamma = 0.0;
alpar@9 978 }
alpar@9 979 else
alpar@9 980 xassert(mip != mip);
alpar@9 981 /* determine the change of x[j] in the adjacent basis:
alpar@9 982 delta x[j] = new x[j] - old x[j] */
alpar@9 983 delta = (kase < 0 ? floor(beta) : ceil(beta)) - beta;
alpar@9 984 /* compute the change of x[k] in the adjacent basis:
alpar@9 985 delta x[k] = new x[k] - old x[k] = delta x[j] / alfa */
alpar@9 986 delta /= alfa;
alpar@9 987 /* compute the change of the objective in the adjacent basis:
alpar@9 988 delta z = new z - old z = gamma * delta x[k] */
alpar@9 989 dz = gamma * delta;
alpar@9 990 if (mip->dir == GLP_MIN)
alpar@9 991 xassert(dz >= 0.0);
alpar@9 992 else if (mip->dir == GLP_MAX)
alpar@9 993 xassert(dz <= 0.0);
alpar@9 994 else
alpar@9 995 xassert(mip != mip);
alpar@9 996 /* compute the new objective value in the adjacent basis:
alpar@9 997 new z = old z + delta z */
alpar@9 998 if (kase < 0)
alpar@9 999 *dn = mip->obj_val + dz;
alpar@9 1000 else
alpar@9 1001 *up = mip->obj_val + dz;
alpar@9 1002 }
alpar@9 1003 /*xprintf("obj = %g; dn = %g; up = %g\n",
alpar@9 1004 mip->obj_val, *dn, *up);*/
alpar@9 1005 return;
alpar@9 1006 }
alpar@9 1007
alpar@9 1008 /***********************************************************************
alpar@9 1009 * NAME
alpar@9 1010 *
alpar@9 1011 * ios_round_bound - improve local bound by rounding
alpar@9 1012 *
alpar@9 1013 * SYNOPSIS
alpar@9 1014 *
alpar@9 1015 * #include "glpios.h"
alpar@9 1016 * double ios_round_bound(glp_tree *tree, double bound);
alpar@9 1017 *
alpar@9 1018 * RETURNS
alpar@9 1019 *
alpar@9 1020 * For the given local bound for any integer feasible solution to the
alpar@9 1021 * current subproblem the routine ios_round_bound returns an improved
alpar@9 1022 * local bound for the same integer feasible solution.
alpar@9 1023 *
alpar@9 1024 * BACKGROUND
alpar@9 1025 *
alpar@9 1026 * Let the current subproblem has the following objective function:
alpar@9 1027 *
alpar@9 1028 * z = sum c[j] * x[j] + s >= b, (1)
alpar@9 1029 * j in J
alpar@9 1030 *
alpar@9 1031 * where J = {j: c[j] is non-zero and integer, x[j] is integer}, s is
alpar@9 1032 * the sum of terms corresponding to fixed variables, b is an initial
alpar@9 1033 * local bound (minimization).
alpar@9 1034 *
alpar@9 1035 * From (1) it follows that:
alpar@9 1036 *
alpar@9 1037 * d * sum (c[j] / d) * x[j] + s >= b, (2)
alpar@9 1038 * j in J
alpar@9 1039 *
alpar@9 1040 * or, equivalently,
alpar@9 1041 *
alpar@9 1042 * sum (c[j] / d) * x[j] >= (b - s) / d = h, (3)
alpar@9 1043 * j in J
alpar@9 1044 *
alpar@9 1045 * where d = gcd(c[j]). Since the left-hand side of (3) is integer,
alpar@9 1046 * h = (b - s) / d can be rounded up to the nearest integer:
alpar@9 1047 *
alpar@9 1048 * h' = ceil(h) = (b' - s) / d, (4)
alpar@9 1049 *
alpar@9 1050 * that gives an rounded, improved local bound:
alpar@9 1051 *
alpar@9 1052 * b' = d * h' + s. (5)
alpar@9 1053 *
alpar@9 1054 * In case of maximization '>=' in (1) should be replaced by '<=' that
alpar@9 1055 * leads to the following formula:
alpar@9 1056 *
alpar@9 1057 * h' = floor(h) = (b' - s) / d, (6)
alpar@9 1058 *
alpar@9 1059 * which should used in the same way as (4).
alpar@9 1060 *
alpar@9 1061 * NOTE: If b is a valid local bound for a child of the current
alpar@9 1062 * subproblem, b' is also valid for that child subproblem. */
alpar@9 1063
alpar@9 1064 double ios_round_bound(glp_tree *tree, double bound)
alpar@9 1065 { glp_prob *mip = tree->mip;
alpar@9 1066 int n = mip->n;
alpar@9 1067 int d, j, nn, *c = tree->iwrk;
alpar@9 1068 double s, h;
alpar@9 1069 /* determine c[j] and compute s */
alpar@9 1070 nn = 0, s = mip->c0, d = 0;
alpar@9 1071 for (j = 1; j <= n; j++)
alpar@9 1072 { GLPCOL *col = mip->col[j];
alpar@9 1073 if (col->coef == 0.0) continue;
alpar@9 1074 if (col->type == GLP_FX)
alpar@9 1075 { /* fixed variable */
alpar@9 1076 s += col->coef * col->prim;
alpar@9 1077 }
alpar@9 1078 else
alpar@9 1079 { /* non-fixed variable */
alpar@9 1080 if (col->kind != GLP_IV) goto skip;
alpar@9 1081 if (col->coef != floor(col->coef)) goto skip;
alpar@9 1082 if (fabs(col->coef) <= (double)INT_MAX)
alpar@9 1083 c[++nn] = (int)fabs(col->coef);
alpar@9 1084 else
alpar@9 1085 d = 1;
alpar@9 1086 }
alpar@9 1087 }
alpar@9 1088 /* compute d = gcd(c[1],...c[nn]) */
alpar@9 1089 if (d == 0)
alpar@9 1090 { if (nn == 0) goto skip;
alpar@9 1091 d = gcdn(nn, c);
alpar@9 1092 }
alpar@9 1093 xassert(d > 0);
alpar@9 1094 /* compute new local bound */
alpar@9 1095 if (mip->dir == GLP_MIN)
alpar@9 1096 { if (bound != +DBL_MAX)
alpar@9 1097 { h = (bound - s) / (double)d;
alpar@9 1098 if (h >= floor(h) + 0.001)
alpar@9 1099 { /* round up */
alpar@9 1100 h = ceil(h);
alpar@9 1101 /*xprintf("d = %d; old = %g; ", d, bound);*/
alpar@9 1102 bound = (double)d * h + s;
alpar@9 1103 /*xprintf("new = %g\n", bound);*/
alpar@9 1104 }
alpar@9 1105 }
alpar@9 1106 }
alpar@9 1107 else if (mip->dir == GLP_MAX)
alpar@9 1108 { if (bound != -DBL_MAX)
alpar@9 1109 { h = (bound - s) / (double)d;
alpar@9 1110 if (h <= ceil(h) - 0.001)
alpar@9 1111 { /* round down */
alpar@9 1112 h = floor(h);
alpar@9 1113 bound = (double)d * h + s;
alpar@9 1114 }
alpar@9 1115 }
alpar@9 1116 }
alpar@9 1117 else
alpar@9 1118 xassert(mip != mip);
alpar@9 1119 skip: return bound;
alpar@9 1120 }
alpar@9 1121
alpar@9 1122 /***********************************************************************
alpar@9 1123 * NAME
alpar@9 1124 *
alpar@9 1125 * ios_is_hopeful - check if subproblem is hopeful
alpar@9 1126 *
alpar@9 1127 * SYNOPSIS
alpar@9 1128 *
alpar@9 1129 * #include "glpios.h"
alpar@9 1130 * int ios_is_hopeful(glp_tree *tree, double bound);
alpar@9 1131 *
alpar@9 1132 * DESCRIPTION
alpar@9 1133 *
alpar@9 1134 * Given the local bound of a subproblem the routine ios_is_hopeful
alpar@9 1135 * checks if the subproblem can have an integer optimal solution which
alpar@9 1136 * is better than the best one currently known.
alpar@9 1137 *
alpar@9 1138 * RETURNS
alpar@9 1139 *
alpar@9 1140 * If the subproblem can have a better integer optimal solution, the
alpar@9 1141 * routine returns non-zero; otherwise, if the corresponding branch can
alpar@9 1142 * be pruned, the routine returns zero. */
alpar@9 1143
alpar@9 1144 int ios_is_hopeful(glp_tree *tree, double bound)
alpar@9 1145 { glp_prob *mip = tree->mip;
alpar@9 1146 int ret = 1;
alpar@9 1147 double eps;
alpar@9 1148 if (mip->mip_stat == GLP_FEAS)
alpar@9 1149 { eps = tree->parm->tol_obj * (1.0 + fabs(mip->mip_obj));
alpar@9 1150 switch (mip->dir)
alpar@9 1151 { case GLP_MIN:
alpar@9 1152 if (bound >= mip->mip_obj - eps) ret = 0;
alpar@9 1153 break;
alpar@9 1154 case GLP_MAX:
alpar@9 1155 if (bound <= mip->mip_obj + eps) ret = 0;
alpar@9 1156 break;
alpar@9 1157 default:
alpar@9 1158 xassert(mip != mip);
alpar@9 1159 }
alpar@9 1160 }
alpar@9 1161 else
alpar@9 1162 { switch (mip->dir)
alpar@9 1163 { case GLP_MIN:
alpar@9 1164 if (bound == +DBL_MAX) ret = 0;
alpar@9 1165 break;
alpar@9 1166 case GLP_MAX:
alpar@9 1167 if (bound == -DBL_MAX) ret = 0;
alpar@9 1168 break;
alpar@9 1169 default:
alpar@9 1170 xassert(mip != mip);
alpar@9 1171 }
alpar@9 1172 }
alpar@9 1173 return ret;
alpar@9 1174 }
alpar@9 1175
alpar@9 1176 /***********************************************************************
alpar@9 1177 * NAME
alpar@9 1178 *
alpar@9 1179 * ios_best_node - find active node with best local bound
alpar@9 1180 *
alpar@9 1181 * SYNOPSIS
alpar@9 1182 *
alpar@9 1183 * #include "glpios.h"
alpar@9 1184 * int ios_best_node(glp_tree *tree);
alpar@9 1185 *
alpar@9 1186 * DESCRIPTION
alpar@9 1187 *
alpar@9 1188 * The routine ios_best_node finds an active node whose local bound is
alpar@9 1189 * best among other active nodes.
alpar@9 1190 *
alpar@9 1191 * It is understood that the integer optimal solution of the original
alpar@9 1192 * mip problem cannot be better than the best bound, so the best bound
alpar@9 1193 * is an lower (minimization) or upper (maximization) global bound for
alpar@9 1194 * the original problem.
alpar@9 1195 *
alpar@9 1196 * RETURNS
alpar@9 1197 *
alpar@9 1198 * The routine ios_best_node returns the subproblem reference number
alpar@9 1199 * for the best node. However, if the tree is empty, it returns zero. */
alpar@9 1200
alpar@9 1201 int ios_best_node(glp_tree *tree)
alpar@9 1202 { IOSNPD *node, *best = NULL;
alpar@9 1203 switch (tree->mip->dir)
alpar@9 1204 { case GLP_MIN:
alpar@9 1205 /* minimization */
alpar@9 1206 for (node = tree->head; node != NULL; node = node->next)
alpar@9 1207 if (best == NULL || best->bound > node->bound)
alpar@9 1208 best = node;
alpar@9 1209 break;
alpar@9 1210 case GLP_MAX:
alpar@9 1211 /* maximization */
alpar@9 1212 for (node = tree->head; node != NULL; node = node->next)
alpar@9 1213 if (best == NULL || best->bound < node->bound)
alpar@9 1214 best = node;
alpar@9 1215 break;
alpar@9 1216 default:
alpar@9 1217 xassert(tree != tree);
alpar@9 1218 }
alpar@9 1219 return best == NULL ? 0 : best->p;
alpar@9 1220 }
alpar@9 1221
alpar@9 1222 /***********************************************************************
alpar@9 1223 * NAME
alpar@9 1224 *
alpar@9 1225 * ios_relative_gap - compute relative mip gap
alpar@9 1226 *
alpar@9 1227 * SYNOPSIS
alpar@9 1228 *
alpar@9 1229 * #include "glpios.h"
alpar@9 1230 * double ios_relative_gap(glp_tree *tree);
alpar@9 1231 *
alpar@9 1232 * DESCRIPTION
alpar@9 1233 *
alpar@9 1234 * The routine ios_relative_gap computes the relative mip gap using the
alpar@9 1235 * formula:
alpar@9 1236 *
alpar@9 1237 * gap = |best_mip - best_bnd| / (|best_mip| + DBL_EPSILON),
alpar@9 1238 *
alpar@9 1239 * where best_mip is the best integer feasible solution found so far,
alpar@9 1240 * best_bnd is the best (global) bound. If no integer feasible solution
alpar@9 1241 * has been found yet, rel_gap is set to DBL_MAX.
alpar@9 1242 *
alpar@9 1243 * RETURNS
alpar@9 1244 *
alpar@9 1245 * The routine ios_relative_gap returns the relative mip gap. */
alpar@9 1246
alpar@9 1247 double ios_relative_gap(glp_tree *tree)
alpar@9 1248 { glp_prob *mip = tree->mip;
alpar@9 1249 int p;
alpar@9 1250 double best_mip, best_bnd, gap;
alpar@9 1251 if (mip->mip_stat == GLP_FEAS)
alpar@9 1252 { best_mip = mip->mip_obj;
alpar@9 1253 p = ios_best_node(tree);
alpar@9 1254 if (p == 0)
alpar@9 1255 { /* the tree is empty */
alpar@9 1256 gap = 0.0;
alpar@9 1257 }
alpar@9 1258 else
alpar@9 1259 { best_bnd = tree->slot[p].node->bound;
alpar@9 1260 gap = fabs(best_mip - best_bnd) / (fabs(best_mip) +
alpar@9 1261 DBL_EPSILON);
alpar@9 1262 }
alpar@9 1263 }
alpar@9 1264 else
alpar@9 1265 { /* no integer feasible solution has been found yet */
alpar@9 1266 gap = DBL_MAX;
alpar@9 1267 }
alpar@9 1268 return gap;
alpar@9 1269 }
alpar@9 1270
alpar@9 1271 /***********************************************************************
alpar@9 1272 * NAME
alpar@9 1273 *
alpar@9 1274 * ios_solve_node - solve LP relaxation of current subproblem
alpar@9 1275 *
alpar@9 1276 * SYNOPSIS
alpar@9 1277 *
alpar@9 1278 * #include "glpios.h"
alpar@9 1279 * int ios_solve_node(glp_tree *tree);
alpar@9 1280 *
alpar@9 1281 * DESCRIPTION
alpar@9 1282 *
alpar@9 1283 * The routine ios_solve_node re-optimizes LP relaxation of the current
alpar@9 1284 * subproblem using the dual simplex method.
alpar@9 1285 *
alpar@9 1286 * RETURNS
alpar@9 1287 *
alpar@9 1288 * The routine returns the code which is reported by glp_simplex. */
alpar@9 1289
alpar@9 1290 int ios_solve_node(glp_tree *tree)
alpar@9 1291 { glp_prob *mip = tree->mip;
alpar@9 1292 glp_smcp parm;
alpar@9 1293 int ret;
alpar@9 1294 /* the current subproblem must exist */
alpar@9 1295 xassert(tree->curr != NULL);
alpar@9 1296 /* set some control parameters */
alpar@9 1297 glp_init_smcp(&parm);
alpar@9 1298 switch (tree->parm->msg_lev)
alpar@9 1299 { case GLP_MSG_OFF:
alpar@9 1300 parm.msg_lev = GLP_MSG_OFF; break;
alpar@9 1301 case GLP_MSG_ERR:
alpar@9 1302 parm.msg_lev = GLP_MSG_ERR; break;
alpar@9 1303 case GLP_MSG_ON:
alpar@9 1304 case GLP_MSG_ALL:
alpar@9 1305 parm.msg_lev = GLP_MSG_ON; break;
alpar@9 1306 case GLP_MSG_DBG:
alpar@9 1307 parm.msg_lev = GLP_MSG_ALL; break;
alpar@9 1308 default:
alpar@9 1309 xassert(tree != tree);
alpar@9 1310 }
alpar@9 1311 parm.meth = GLP_DUALP;
alpar@9 1312 if (tree->parm->msg_lev < GLP_MSG_DBG)
alpar@9 1313 parm.out_dly = tree->parm->out_dly;
alpar@9 1314 else
alpar@9 1315 parm.out_dly = 0;
alpar@9 1316 /* if the incumbent objective value is already known, use it to
alpar@9 1317 prematurely terminate the dual simplex search */
alpar@9 1318 if (mip->mip_stat == GLP_FEAS)
alpar@9 1319 { switch (tree->mip->dir)
alpar@9 1320 { case GLP_MIN:
alpar@9 1321 parm.obj_ul = mip->mip_obj;
alpar@9 1322 break;
alpar@9 1323 case GLP_MAX:
alpar@9 1324 parm.obj_ll = mip->mip_obj;
alpar@9 1325 break;
alpar@9 1326 default:
alpar@9 1327 xassert(mip != mip);
alpar@9 1328 }
alpar@9 1329 }
alpar@9 1330 /* try to solve/re-optimize the LP relaxation */
alpar@9 1331 ret = glp_simplex(mip, &parm);
alpar@9 1332 tree->curr->solved++;
alpar@9 1333 #if 0
alpar@9 1334 xprintf("ret = %d; status = %d; pbs = %d; dbs = %d; some = %d\n",
alpar@9 1335 ret, glp_get_status(mip), mip->pbs_stat, mip->dbs_stat,
alpar@9 1336 mip->some);
alpar@9 1337 lpx_print_sol(mip, "sol");
alpar@9 1338 #endif
alpar@9 1339 return ret;
alpar@9 1340 }
alpar@9 1341
alpar@9 1342 /**********************************************************************/
alpar@9 1343
alpar@9 1344 IOSPOOL *ios_create_pool(glp_tree *tree)
alpar@9 1345 { /* create cut pool */
alpar@9 1346 IOSPOOL *pool;
alpar@9 1347 #if 0
alpar@9 1348 pool = dmp_get_atom(tree->pool, sizeof(IOSPOOL));
alpar@9 1349 #else
alpar@9 1350 xassert(tree == tree);
alpar@9 1351 pool = xmalloc(sizeof(IOSPOOL));
alpar@9 1352 #endif
alpar@9 1353 pool->size = 0;
alpar@9 1354 pool->head = pool->tail = NULL;
alpar@9 1355 pool->ord = 0, pool->curr = NULL;
alpar@9 1356 return pool;
alpar@9 1357 }
alpar@9 1358
alpar@9 1359 int ios_add_row(glp_tree *tree, IOSPOOL *pool,
alpar@9 1360 const char *name, int klass, int flags, int len, const int ind[],
alpar@9 1361 const double val[], int type, double rhs)
alpar@9 1362 { /* add row (constraint) to the cut pool */
alpar@9 1363 IOSCUT *cut;
alpar@9 1364 IOSAIJ *aij;
alpar@9 1365 int k;
alpar@9 1366 xassert(pool != NULL);
alpar@9 1367 cut = dmp_get_atom(tree->pool, sizeof(IOSCUT));
alpar@9 1368 if (name == NULL || name[0] == '\0')
alpar@9 1369 cut->name = NULL;
alpar@9 1370 else
alpar@9 1371 { for (k = 0; name[k] != '\0'; k++)
alpar@9 1372 { if (k == 256)
alpar@9 1373 xerror("glp_ios_add_row: cut name too long\n");
alpar@9 1374 if (iscntrl((unsigned char)name[k]))
alpar@9 1375 xerror("glp_ios_add_row: cut name contains invalid chara"
alpar@9 1376 "cter(s)\n");
alpar@9 1377 }
alpar@9 1378 cut->name = dmp_get_atom(tree->pool, strlen(name)+1);
alpar@9 1379 strcpy(cut->name, name);
alpar@9 1380 }
alpar@9 1381 if (!(0 <= klass && klass <= 255))
alpar@9 1382 xerror("glp_ios_add_row: klass = %d; invalid cut class\n",
alpar@9 1383 klass);
alpar@9 1384 cut->klass = (unsigned char)klass;
alpar@9 1385 if (flags != 0)
alpar@9 1386 xerror("glp_ios_add_row: flags = %d; invalid cut flags\n",
alpar@9 1387 flags);
alpar@9 1388 cut->ptr = NULL;
alpar@9 1389 if (!(0 <= len && len <= tree->n))
alpar@9 1390 xerror("glp_ios_add_row: len = %d; invalid cut length\n",
alpar@9 1391 len);
alpar@9 1392 for (k = 1; k <= len; k++)
alpar@9 1393 { aij = dmp_get_atom(tree->pool, sizeof(IOSAIJ));
alpar@9 1394 if (!(1 <= ind[k] && ind[k] <= tree->n))
alpar@9 1395 xerror("glp_ios_add_row: ind[%d] = %d; column index out of "
alpar@9 1396 "range\n", k, ind[k]);
alpar@9 1397 aij->j = ind[k];
alpar@9 1398 aij->val = val[k];
alpar@9 1399 aij->next = cut->ptr;
alpar@9 1400 cut->ptr = aij;
alpar@9 1401 }
alpar@9 1402 if (!(type == GLP_LO || type == GLP_UP || type == GLP_FX))
alpar@9 1403 xerror("glp_ios_add_row: type = %d; invalid cut type\n",
alpar@9 1404 type);
alpar@9 1405 cut->type = (unsigned char)type;
alpar@9 1406 cut->rhs = rhs;
alpar@9 1407 cut->prev = pool->tail;
alpar@9 1408 cut->next = NULL;
alpar@9 1409 if (cut->prev == NULL)
alpar@9 1410 pool->head = cut;
alpar@9 1411 else
alpar@9 1412 cut->prev->next = cut;
alpar@9 1413 pool->tail = cut;
alpar@9 1414 pool->size++;
alpar@9 1415 return pool->size;
alpar@9 1416 }
alpar@9 1417
alpar@9 1418 IOSCUT *ios_find_row(IOSPOOL *pool, int i)
alpar@9 1419 { /* find row (constraint) in the cut pool */
alpar@9 1420 /* (smart linear search) */
alpar@9 1421 xassert(pool != NULL);
alpar@9 1422 xassert(1 <= i && i <= pool->size);
alpar@9 1423 if (pool->ord == 0)
alpar@9 1424 { xassert(pool->curr == NULL);
alpar@9 1425 pool->ord = 1;
alpar@9 1426 pool->curr = pool->head;
alpar@9 1427 }
alpar@9 1428 xassert(pool->curr != NULL);
alpar@9 1429 if (i < pool->ord)
alpar@9 1430 { if (i < pool->ord - i)
alpar@9 1431 { pool->ord = 1;
alpar@9 1432 pool->curr = pool->head;
alpar@9 1433 while (pool->ord != i)
alpar@9 1434 { pool->ord++;
alpar@9 1435 xassert(pool->curr != NULL);
alpar@9 1436 pool->curr = pool->curr->next;
alpar@9 1437 }
alpar@9 1438 }
alpar@9 1439 else
alpar@9 1440 { while (pool->ord != i)
alpar@9 1441 { pool->ord--;
alpar@9 1442 xassert(pool->curr != NULL);
alpar@9 1443 pool->curr = pool->curr->prev;
alpar@9 1444 }
alpar@9 1445 }
alpar@9 1446 }
alpar@9 1447 else if (i > pool->ord)
alpar@9 1448 { if (i - pool->ord < pool->size - i)
alpar@9 1449 { while (pool->ord != i)
alpar@9 1450 { pool->ord++;
alpar@9 1451 xassert(pool->curr != NULL);
alpar@9 1452 pool->curr = pool->curr->next;
alpar@9 1453 }
alpar@9 1454 }
alpar@9 1455 else
alpar@9 1456 { pool->ord = pool->size;
alpar@9 1457 pool->curr = pool->tail;
alpar@9 1458 while (pool->ord != i)
alpar@9 1459 { pool->ord--;
alpar@9 1460 xassert(pool->curr != NULL);
alpar@9 1461 pool->curr = pool->curr->prev;
alpar@9 1462 }
alpar@9 1463 }
alpar@9 1464 }
alpar@9 1465 xassert(pool->ord == i);
alpar@9 1466 xassert(pool->curr != NULL);
alpar@9 1467 return pool->curr;
alpar@9 1468 }
alpar@9 1469
alpar@9 1470 void ios_del_row(glp_tree *tree, IOSPOOL *pool, int i)
alpar@9 1471 { /* remove row (constraint) from the cut pool */
alpar@9 1472 IOSCUT *cut;
alpar@9 1473 IOSAIJ *aij;
alpar@9 1474 xassert(pool != NULL);
alpar@9 1475 if (!(1 <= i && i <= pool->size))
alpar@9 1476 xerror("glp_ios_del_row: i = %d; cut number out of range\n",
alpar@9 1477 i);
alpar@9 1478 cut = ios_find_row(pool, i);
alpar@9 1479 xassert(pool->curr == cut);
alpar@9 1480 if (cut->next != NULL)
alpar@9 1481 pool->curr = cut->next;
alpar@9 1482 else if (cut->prev != NULL)
alpar@9 1483 pool->ord--, pool->curr = cut->prev;
alpar@9 1484 else
alpar@9 1485 pool->ord = 0, pool->curr = NULL;
alpar@9 1486 if (cut->name != NULL)
alpar@9 1487 dmp_free_atom(tree->pool, cut->name, strlen(cut->name)+1);
alpar@9 1488 if (cut->prev == NULL)
alpar@9 1489 { xassert(pool->head == cut);
alpar@9 1490 pool->head = cut->next;
alpar@9 1491 }
alpar@9 1492 else
alpar@9 1493 { xassert(cut->prev->next == cut);
alpar@9 1494 cut->prev->next = cut->next;
alpar@9 1495 }
alpar@9 1496 if (cut->next == NULL)
alpar@9 1497 { xassert(pool->tail == cut);
alpar@9 1498 pool->tail = cut->prev;
alpar@9 1499 }
alpar@9 1500 else
alpar@9 1501 { xassert(cut->next->prev == cut);
alpar@9 1502 cut->next->prev = cut->prev;
alpar@9 1503 }
alpar@9 1504 while (cut->ptr != NULL)
alpar@9 1505 { aij = cut->ptr;
alpar@9 1506 cut->ptr = aij->next;
alpar@9 1507 dmp_free_atom(tree->pool, aij, sizeof(IOSAIJ));
alpar@9 1508 }
alpar@9 1509 dmp_free_atom(tree->pool, cut, sizeof(IOSCUT));
alpar@9 1510 pool->size--;
alpar@9 1511 return;
alpar@9 1512 }
alpar@9 1513
alpar@9 1514 void ios_clear_pool(glp_tree *tree, IOSPOOL *pool)
alpar@9 1515 { /* remove all rows (constraints) from the cut pool */
alpar@9 1516 xassert(pool != NULL);
alpar@9 1517 while (pool->head != NULL)
alpar@9 1518 { IOSCUT *cut = pool->head;
alpar@9 1519 pool->head = cut->next;
alpar@9 1520 if (cut->name != NULL)
alpar@9 1521 dmp_free_atom(tree->pool, cut->name, strlen(cut->name)+1);
alpar@9 1522 while (cut->ptr != NULL)
alpar@9 1523 { IOSAIJ *aij = cut->ptr;
alpar@9 1524 cut->ptr = aij->next;
alpar@9 1525 dmp_free_atom(tree->pool, aij, sizeof(IOSAIJ));
alpar@9 1526 }
alpar@9 1527 dmp_free_atom(tree->pool, cut, sizeof(IOSCUT));
alpar@9 1528 }
alpar@9 1529 pool->size = 0;
alpar@9 1530 pool->head = pool->tail = NULL;
alpar@9 1531 pool->ord = 0, pool->curr = NULL;
alpar@9 1532 return;
alpar@9 1533 }
alpar@9 1534
alpar@9 1535 void ios_delete_pool(glp_tree *tree, IOSPOOL *pool)
alpar@9 1536 { /* delete cut pool */
alpar@9 1537 xassert(pool != NULL);
alpar@9 1538 ios_clear_pool(tree, pool);
alpar@9 1539 xfree(pool);
alpar@9 1540 return;
alpar@9 1541 }
alpar@9 1542
alpar@9 1543 /**********************************************************************/
alpar@9 1544
alpar@9 1545 #if 0
alpar@9 1546 static int refer_to_node(glp_tree *tree, int j)
alpar@9 1547 { /* determine node number corresponding to binary variable x[j] or
alpar@9 1548 its complement */
alpar@9 1549 glp_prob *mip = tree->mip;
alpar@9 1550 int n = mip->n;
alpar@9 1551 int *ref;
alpar@9 1552 if (j > 0)
alpar@9 1553 ref = tree->n_ref;
alpar@9 1554 else
alpar@9 1555 ref = tree->c_ref, j = - j;
alpar@9 1556 xassert(1 <= j && j <= n);
alpar@9 1557 if (ref[j] == 0)
alpar@9 1558 { /* new node is needed */
alpar@9 1559 SCG *g = tree->g;
alpar@9 1560 int n_max = g->n_max;
alpar@9 1561 ref[j] = scg_add_nodes(g, 1);
alpar@9 1562 if (g->n_max > n_max)
alpar@9 1563 { int *save = tree->j_ref;
alpar@9 1564 tree->j_ref = xcalloc(1+g->n_max, sizeof(int));
alpar@9 1565 memcpy(&tree->j_ref[1], &save[1], g->n * sizeof(int));
alpar@9 1566 xfree(save);
alpar@9 1567 }
alpar@9 1568 xassert(ref[j] == g->n);
alpar@9 1569 tree->j_ref[ref[j]] = j;
alpar@9 1570 xassert(tree->curr != NULL);
alpar@9 1571 if (tree->curr->level > 0) tree->curr->own_nn++;
alpar@9 1572 }
alpar@9 1573 return ref[j];
alpar@9 1574 }
alpar@9 1575 #endif
alpar@9 1576
alpar@9 1577 #if 0
alpar@9 1578 void ios_add_edge(glp_tree *tree, int j1, int j2)
alpar@9 1579 { /* add new edge to the conflict graph */
alpar@9 1580 glp_prob *mip = tree->mip;
alpar@9 1581 int n = mip->n;
alpar@9 1582 SCGRIB *e;
alpar@9 1583 int first, i1, i2;
alpar@9 1584 xassert(-n <= j1 && j1 <= +n && j1 != 0);
alpar@9 1585 xassert(-n <= j2 && j2 <= +n && j2 != 0);
alpar@9 1586 xassert(j1 != j2);
alpar@9 1587 /* determine number of the first node, which was added for the
alpar@9 1588 current subproblem */
alpar@9 1589 xassert(tree->curr != NULL);
alpar@9 1590 first = tree->g->n - tree->curr->own_nn + 1;
alpar@9 1591 /* determine node numbers for both endpoints */
alpar@9 1592 i1 = refer_to_node(tree, j1);
alpar@9 1593 i2 = refer_to_node(tree, j2);
alpar@9 1594 /* add edge (i1,i2) to the conflict graph */
alpar@9 1595 e = scg_add_edge(tree->g, i1, i2);
alpar@9 1596 /* if the current subproblem is not the root and both endpoints
alpar@9 1597 were created on some previous levels, save the edge */
alpar@9 1598 if (tree->curr->level > 0 && i1 < first && i2 < first)
alpar@9 1599 { IOSRIB *rib;
alpar@9 1600 rib = dmp_get_atom(tree->pool, sizeof(IOSRIB));
alpar@9 1601 rib->j1 = j1;
alpar@9 1602 rib->j2 = j2;
alpar@9 1603 rib->e = e;
alpar@9 1604 rib->next = tree->curr->e_ptr;
alpar@9 1605 tree->curr->e_ptr = rib;
alpar@9 1606 }
alpar@9 1607 return;
alpar@9 1608 }
alpar@9 1609 #endif
alpar@9 1610
alpar@9 1611 /* eof */