3
/*----------------------------------------------------------------------
4
-- This code is part of GNU Linear Programming Kit (GLPK).
6
-- Copyright (C) 2000, 01, 02, 03, 04, 05, 06 Andrew Makhorin,
7
-- Department for Applied Informatics, Moscow Aviation Institute,
8
-- Moscow, Russia. All rights reserved. E-mail: <mao@mai2.rcnet.ru>.
10
-- GLPK is free software; you can redistribute it and/or modify it
11
-- under the terms of the GNU General Public License as published by
12
-- the Free Software Foundation; either version 2, or (at your option)
15
-- GLPK is distributed in the hope that it will be useful, but WITHOUT
16
-- ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
17
-- or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
18
-- License for more details.
20
-- You should have received a copy of the GNU General Public License
21
-- along with GLPK; see the file COPYING. If not, write to the Free
22
-- Software Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA
24
----------------------------------------------------------------------*/
32
/*----------------------------------------------------------------------
33
-- ipp_create_wksp - create MIP presolver workspace.
37
-- #include "glpipp.h"
38
-- IPP *ipp_create_wksp(void);
42
-- The routine ipp_create_wksp creates an empty workspace used by the
43
-- MIP presolver routines.
47
-- The routine returns a pointer to the MIP workspace created. */
49
IPP *ipp_create_wksp(void)
51
ipp = umalloc(sizeof(IPP));
55
ipp->orig_dir = LPX_MIN;
57
ipp->row_pool = dmp_create_pool(sizeof(IPPROW));
58
ipp->col_pool = dmp_create_pool(sizeof(IPPCOL));
59
ipp->aij_pool = dmp_create_pool(sizeof(IPPAIJ));
65
ipp->tqe_pool = dmp_create_pool(0);
72
/*----------------------------------------------------------------------
73
-- ipp_add_row - add new row to the transformed problem.
77
-- #include "glpipp.h"
78
-- IPPROW *ipp_add_row(IPP *ipp, gnm_float lb, gnm_float ub);
82
-- The routine ipp_add_row adds a new empty row to the transformed
85
-- The parameter lb is an lower bound of the new row (-DBL_MAX means
86
-- the row has no lower bound).
88
-- The parameter ub is an upper bound of the new row (+DBL_MAX means
89
-- the row has no upper bound).
93
-- The routine returns a pointer to the created row. */
95
IPPROW *ipp_add_row(IPP *ipp, gnm_float lb, gnm_float ub)
97
/* perform sanity checks */
100
row = dmp_get_atom(ipp->row_pool);
106
row->next = ipp->row_ptr;
110
/* add the row to the linked list of rows */
111
if (row->next != NULL) row->next->prev = row;
116
/*----------------------------------------------------------------------
117
-- ipp_add_col - add new column to the transformed problem.
121
-- #include "glpipp.h"
122
-- IPPCOL *ipp_add_col(IPP *ipp, int i_flag, gnm_float lb, gnm_float ub,
127
-- The routine ipp_add_col adds a new empty column to the transformed
130
-- The parameter i_flag is a column integrality flag. If it is set, the
131
-- new column is integer, otherwise the new column is continuous.
133
-- The parameter lb is an lower bound of the new column (-DBL_MAX means
134
-- the column has no lower bound).
136
-- The parameter ub is an upper bound of the new column (+DBL_MAX means
137
-- the column has no upper bound).
139
-- The parameter c is an objective coefficient at the new column.
143
-- The routine returns a pointer to the created column. */
145
IPPCOL *ipp_add_col(IPP *ipp, int i_flag, gnm_float lb, gnm_float ub,
148
/* perform sanity checks */
151
{ if (lb != -DBL_MAX) insist(lb == gnm_floor(lb));
152
if (ub != +DBL_MAX) insist(ub == gnm_floor(ub));
154
/* create new column */
155
col = dmp_get_atom(ipp->col_pool);
156
col->j = ++(ipp->ncols);
157
col->i_flag = i_flag;
164
col->next = ipp->col_ptr;
168
/* add the column to the linked list of columns */
169
if (col->next != NULL) col->next->prev = col;
174
/*----------------------------------------------------------------------
175
-- ipp_add_aij - add new element to the constraint matrix.
179
-- #include "glpipp.h"
180
-- IPPAIJ *ipp_add_aij(IPP *ipp, IPPROW *row, IPPCOL *col, gnm_float val);
184
-- The routine ipp_add_aij adds a new element to the constraint matrix
185
-- of the transformed problem.
187
-- The parameter row is a pointer to a row, in which the new element
190
-- The parameter col is a pointer to a column, in which the new element
193
-- The parameter val is a numeric value of the new element.
197
-- The routine returns a pointer to the created element. */
199
IPPAIJ *ipp_add_aij(IPP *ipp, IPPROW *row, IPPCOL *col, gnm_float val)
202
aij = dmp_get_atom(ipp->aij_pool);
207
aij->r_next = row->ptr;
209
aij->c_next = col->ptr;
210
if (row->ptr != NULL) row->ptr->r_prev = aij;
211
if (col->ptr != NULL) col->ptr->c_prev = aij;
212
row->ptr = col->ptr = aij;
216
/*----------------------------------------------------------------------
217
-- ipp_remove_row - remove row from the transformed problem.
221
-- #include "glpipp.h"
222
-- void ipp_remove_row(IPP *ipp, IPPROW *row);
226
-- The routine ipp_remove_row removes a row, which the parameter row
227
-- points to, from the transformed problem. */
229
void ipp_remove_row(IPP *ipp, IPPROW *row)
231
/* remove the row from the active queue */
232
ipp_deque_row(ipp, row);
233
/* remove elements of the row from the constraint matrix */
234
while (row->ptr != NULL)
235
{ /* get a next element in the row */
237
/* remove the element from the row list */
238
row->ptr = aij->r_next;
239
/* remove the element from the column list */
240
if (aij->c_prev == NULL)
241
aij->col->ptr = aij->c_next;
243
aij->c_prev->c_next = aij->c_next;
244
if (aij->c_next == NULL)
247
aij->c_next->c_prev = aij->c_prev;
248
/* and return it to its pool */
249
dmp_free_atom(ipp->aij_pool, aij);
251
/* remove the row from the linked list */
252
if (row->prev == NULL)
253
ipp->row_ptr = row->next;
255
row->prev->next = row->next;
256
if (row->next == NULL)
259
row->next->prev = row->prev;
260
/* and return the row to its pool */
261
dmp_free_atom(ipp->row_pool, row);
265
/*----------------------------------------------------------------------
266
-- ipp_remove_col - remove column from the transformed problem.
270
-- #include "glpipp.h"
271
-- void ipp_remove_col(IPP *ipp, IPPCOL *col);
275
-- The routine ipp_remove_col removes a column, which the parameter col
276
-- points to, from the transformed problem. */
278
void ipp_remove_col(IPP *ipp, IPPCOL *col)
280
/* remove the column from the active queue */
281
ipp_deque_col(ipp, col);
282
/* remove elements of the column from the constraint matrix */
283
while (col->ptr != NULL)
284
{ /* get a next element in the column */
286
/* remove the element from the column list */
287
col->ptr = aij->c_next;
288
/* remove the element from the row list */
289
if (aij->r_prev == NULL)
290
aij->row->ptr = aij->r_next;
292
aij->r_prev->r_next = aij->r_next;
293
if (aij->r_next == NULL)
296
aij->r_next->r_prev = aij->r_prev;
297
/* and return it to its pool */
298
dmp_free_atom(ipp->aij_pool, aij);
300
/* remove the column from the linked list */
301
if (col->prev == NULL)
302
ipp->col_ptr = col->next;
304
col->prev->next = col->next;
305
if (col->next == NULL)
308
col->next->prev = col->prev;
309
/* and return the column to its pool */
310
dmp_free_atom(ipp->col_pool, col);
314
/*----------------------------------------------------------------------
315
-- ipp_enque_row - place row in the active queue.
319
-- #include "glpipp.h"
320
-- void ipp_enque_row(IPP *ipp, IPPROW *row);
324
-- The routine ipp_enque_row places the specified row to the queue of
327
void ipp_enque_row(IPP *ipp, IPPROW *row)
331
row->q_next = ipp->row_que;
332
if (ipp->row_que != NULL) ipp->row_que->q_prev = row;
338
/*----------------------------------------------------------------------
339
-- ipp_deque_row - remove row from the active queue.
343
-- #include "glpipp.h"
344
-- void ipp_deque_row(IPP *ipp, IPPROW *row);
348
-- The routine ipp_deque_row removes the specified row from the queue
349
-- of active rows. */
351
void ipp_deque_row(IPP *ipp, IPPROW *row)
354
if (row->q_prev == NULL)
355
ipp->row_que = row->q_next;
357
row->q_prev->q_next = row->q_next;
358
if (row->q_next == NULL)
361
row->q_next->q_prev = row->q_prev;
366
/*----------------------------------------------------------------------
367
-- ipp_enque_col - place column in the active queue.
371
-- #include "glpipp.h"
372
-- void ipp_enque_col(IPP *ipp, IPPCOL *col);
376
-- The routine ipp_enque_col places the specified column to the queue of
377
-- active columns. */
379
void ipp_enque_col(IPP *ipp, IPPCOL *col)
383
col->q_next = ipp->col_que;
384
if (ipp->col_que != NULL) ipp->col_que->q_prev = col;
390
/*----------------------------------------------------------------------
391
-- ipp_deque_col - remove column from the active queue.
395
-- #include "glpipp.h"
396
-- void ipp_deque_col(IPP *ipp, IPPCOL *col);
400
-- The routine ipp_deque_col removes the specified column from the queue
401
-- of active columns. */
403
void ipp_deque_col(IPP *ipp, IPPCOL *col)
406
if (col->q_prev == NULL)
407
ipp->col_que = col->q_next;
409
col->q_prev->q_next = col->q_next;
410
if (col->q_next == NULL)
413
col->q_next->q_prev = col->q_prev;
418
/*----------------------------------------------------------------------
419
-- ipp_append_tqe - append new transformation queue entry.
423
-- #include "glpipp.h"
424
-- void *ipp_append_tqe(IPP *ipp, int type, int size);
428
-- The routine ipp_append_tqe appends a new transformation queue entry
431
-- The parameter type is the entry type.
433
-- The parameter size is the size of a specific part, in bytes.
437
-- The routine returns a pointer to a specific part, which is allocated
438
-- and attached to the entry. */
440
void *ipp_append_tqe(IPP *ipp, int type, int size)
442
tqe = dmp_get_atomv(ipp->tqe_pool, sizeof(IPPTQE));
444
tqe->info = dmp_get_atomv(ipp->tqe_pool, size);
445
tqe->next = ipp->tqe_list;
450
/*----------------------------------------------------------------------
451
-- ipp_load_orig - load original problem into MIP presolver workspace.
455
-- #include "glpipp.h"
456
-- void ipp_load_orig(IPP *ipp, LPX *orig);
460
-- The routine ipp_load_orig loads an original MIP problem, which the
461
-- parameter orig points to, into the MIP presolver workspace.
463
-- On exit from the routine the tranformed problem in the workspace is
464
-- identical to the original problem. */
466
void ipp_load_orig(IPP *ipp, LPX *orig)
469
int i, j, k, type, len, *ind;
470
gnm_float lb, ub, *val;
471
/* save some information about the original problem */
472
ipp->orig_m = lpx_get_num_rows(orig);
473
ipp->orig_n = lpx_get_num_cols(orig);
474
ipp->orig_nnz = lpx_get_num_nz(orig);
475
ipp->orig_dir = lpx_get_obj_dir(orig);
476
/* allocate working arrays */
477
row = ucalloc(1+ipp->orig_m, sizeof(IPPROW *));
478
ind = ucalloc(1+ipp->orig_m, sizeof(int));
479
val = ucalloc(1+ipp->orig_m, sizeof(gnm_float));
480
/* copy rows of the original problem into the workspace */
481
for (i = 1; i <= ipp->orig_m; i++)
482
{ type = lpx_get_row_type(orig, i);
483
if (type == LPX_FR || type == LPX_UP)
486
lb = lpx_get_row_lb(orig, i);
487
if (type == LPX_FR || type == LPX_LO)
490
ub = lpx_get_row_ub(orig, i);
491
row[i] = ipp_add_row(ipp, lb, ub);
493
/* copy columns of the original problem into the workspace; each
494
column created in the workspace is assigned a reference number
495
which is its ordinal number in the original problem */
496
for (j = 1; j <= ipp->orig_n; j++)
497
{ type = lpx_get_col_type(orig, j);
498
if (type == LPX_FR || type == LPX_UP)
501
lb = lpx_get_col_lb(orig, j);
502
if (type == LPX_FR || type == LPX_LO)
505
ub = lpx_get_col_ub(orig, j);
506
col = ipp_add_col(ipp, lpx_get_col_kind(orig, j) == LPX_IV,
507
lb, ub, lpx_get_obj_coef(orig, j));
508
len = lpx_get_mat_col(orig, j, ind, val);
509
for (k = 1; k <= len; k++)
510
ipp_add_aij(ipp, row[ind[k]], col, val[k]);
512
/* copy the constant term of the original objective function */
513
ipp->c0 = lpx_get_obj_coef(orig, 0);
514
/* if the original problem is maximization, change the sign of
515
the objective function, because the transformed problem to be
516
processed by the presolver must be minimization */
517
if (ipp->orig_dir == LPX_MAX)
518
{ for (col = ipp->col_ptr; col != NULL; col = col->next)
522
/* free working arrays */
529
/*----------------------------------------------------------------------
530
-- ipp_tight_bnds - tight current column bounds using implied bounds.
534
-- #include "glpipp.h"
535
-- int ipp_tight_bnds(IPP *ipp, IPPCOL *col, gnm_float lb, gnm_float ub);
539
-- The routines ipp_tight_bnds replaces current bounds of the column,
540
-- which the parameter col points to:
546
-- max(l, l') <= x <= min(u, u'), (2)
548
-- where l' and u' are specified implied bounds of the column.
552
-- 0 - bounds remain unchanged
553
-- 1 - bounds have been changed
554
-- 2 - new bounds are primal infeasible */
556
int ipp_tight_bnds(IPP *ipp, IPPCOL *col, gnm_float lb, gnm_float ub)
560
/* if the column is integral, round implied bounds */
562
{ eps = 1e-5 * (1.0 + gnm_abs(lb));
563
if (gnm_abs(lb - gnm_floor(lb + 0.5)) <= eps)
564
lb = gnm_floor(lb + 0.5);
567
eps = 1e-5 * (1.0 + gnm_abs(ub));
568
if (gnm_abs(ub - gnm_floor(ub + 0.5)) <= eps)
569
ub = gnm_floor(ub + 0.5);
573
/* check for primal infeasibility */
574
if (col->lb != -DBL_MAX)
575
{ eps = 1e-5 * (1.0 + gnm_abs(col->lb));
576
if (ub < col->lb - eps)
581
if (col->ub != +DBL_MAX)
582
{ eps = 1e-5 * (1.0 + gnm_abs(col->ub));
583
if (lb > col->ub + eps)
589
if (col->i_flag && lb > ub + 0.5)
590
{ /* this may happen due to rounding implied bounds */
595
/* replace current bounds by implied ones, if necessary */
597
{ eps = 1e-7 * (1.0 + gnm_abs(lb));
598
if (col->lb < lb - eps)
604
{ eps = 1e-7 * (1.0 + gnm_abs(ub));
605
if (col->ub > ub + eps)
610
/* if new bounds of the column are close to each other, the
611
column can be fixed */
612
if (ret == 1 && col->lb != -DBL_MAX && col->ub != +DBL_MAX)
613
{ eps = 1e-7 * (1.0 + gnm_abs(col->lb));
614
if (col->lb >= col->ub - eps)
615
{ if (gnm_abs(col->lb) <= gnm_abs(col->ub))
621
done: insist(col->lb <= col->ub);
625
/*----------------------------------------------------------------------
626
-- ipp_build_prob - build resultant problem.
630
-- #include "glpipp.h"
631
-- LPX *ipp_build_prob(IPP *ipp);
635
-- The routine ipp_build_prob converts the resultant MIP problem from
636
-- an internal format, in which the problem is stored in the workspace,
637
-- to the standard problem object.
641
-- The routine returns a pointer to the problem object. */
643
LPX *ipp_build_prob(IPP *ipp)
648
int i, j, type, len, *ind;
650
/* create problem object */
651
prob = lpx_create_prob();
652
lpx_set_class(prob, LPX_MIP);
653
/* the resultant problem should have the same optimization sense
654
as the original problem */
655
lpx_set_obj_dir(prob, ipp->orig_dir);
656
/* set the constant term of the objective function */
657
lpx_set_obj_coef(prob, 0,
658
ipp->orig_dir == LPX_MIN ? + ipp->c0 : - ipp->c0);
659
/* copy rows of the resultant problem */
660
for (row = ipp->row_ptr; row != NULL; row = row->next)
661
{ i = lpx_add_rows(prob, 1);
662
if (row->lb == -DBL_MAX && row->ub == +DBL_MAX)
664
else if (row->ub == +DBL_MAX)
666
else if (row->lb == -DBL_MAX)
668
else if (row->lb != row->ub)
672
lpx_set_row_bnds(prob, i, type, row->lb, row->ub);
675
/* copy columns of the resultant problem */
676
ind = ucalloc(1+lpx_get_num_rows(prob), sizeof(int));
677
val = ucalloc(1+lpx_get_num_rows(prob), sizeof(gnm_float));
678
for (col = ipp->col_ptr; col != NULL; col = col->next)
679
{ j = lpx_add_cols(prob, 1);
680
if (col->i_flag) lpx_set_col_kind(prob, j, LPX_IV);
681
if (col->lb == -DBL_MAX && col->ub == +DBL_MAX)
683
else if (col->ub == +DBL_MAX)
685
else if (col->lb == -DBL_MAX)
687
else if (col->lb != col->ub)
691
lpx_set_col_bnds(prob, j, type, col->lb, col->ub);
692
lpx_set_obj_coef(prob, j,
693
ipp->orig_dir == LPX_MIN ? + col->c : - col->c);
694
/* copy constraint coefficients */
696
for (aij = col->ptr; aij != NULL; aij = aij->c_next)
698
ind[len] = aij->row->temp;
701
lpx_set_mat_col(prob, j, len, ind, val);
708
/*----------------------------------------------------------------------
709
-- ipp_load_sol - load solution into MIP presolver workspace.
713
-- #include "glpipp.h"
714
-- void ipp_load_sol(IPP *ipp, LPX *prob);
718
-- The routine ipp_load_sol loads an integer solution of the resultant
719
-- MIP problem into the MIP presolver workspace. */
721
void ipp_load_sol(IPP *ipp, LPX *prob)
724
insist(lpx_mip_status(prob) != LPX_I_UNDEF);
725
ipp->col_stat = ucalloc(1+ipp->ncols, sizeof(int));
726
ipp->col_mipx = ucalloc(1+ipp->ncols, sizeof(gnm_float));
727
for (j = 1; j <= ipp->ncols; j++) ipp->col_stat[j] = 0;
728
/* columns in the problem object follow in the same order as in
729
the column list (see ipp_build_prob) */
731
for (col = ipp->col_ptr; col != NULL; col = col->next)
733
ipp->col_stat[col->j] = 1;
734
ipp->col_mipx[col->j] = lpx_mip_col_val(prob, j);
739
/*----------------------------------------------------------------------
740
-- ipp_unload_sol - unload solution from MIP presolver workspace.
744
-- #include "glpipp.h"
745
-- void ipp_unload_sol(IPP *ipp, LPX *orig, int i_stat);
749
-- The routine ipp_unload_sol unloads a recovered solution from the
750
-- MIP presolver workspace into the original problem object, which the
751
-- parameter orig points to. */
753
void ipp_unload_sol(IPP *ipp, LPX *orig, int i_stat)
754
{ int i, j, k, len, *ind;
755
gnm_float temp, *row_mipx, *val;
756
insist(ipp->orig_m == lpx_get_num_rows(orig));
757
insist(ipp->orig_n == lpx_get_num_cols(orig));
758
insist(ipp->orig_dir == lpx_get_obj_dir(orig));
759
/* all columns must be computed/recovered */
760
insist(ipp->orig_n <= ipp->ncols);
761
for (j = 1; j <= ipp->ncols; j++) insist(ipp->col_stat[j]);
762
/* compute values of auxiliary variables using known values of
763
structural variables (columns) */
764
row_mipx = ucalloc(1+ipp->orig_m, sizeof(gnm_float));
765
ind = ucalloc(1+ipp->orig_n, sizeof(int));
766
val = ucalloc(1+ipp->orig_n, sizeof(gnm_float));
767
for (i = 1; i <= ipp->orig_m; i++)
768
{ len = lpx_get_mat_row(orig, i, ind, val);
770
for (k = 1; k <= len; k++)
771
temp += val[k] * ipp->col_mipx[ind[k]];
776
/* store solution components into the original problem object */
777
lpx_put_mip_soln(orig, i_stat, row_mipx, ipp->col_mipx);
782
/*----------------------------------------------------------------------
783
-- ipp_delete_wksp - delete MIP presolver workspace.
787
-- #include "glpipp.h"
788
-- void ipp_delete_wksp(IPP *ipp);
792
-- The routine ipp_delete_wksp deletes a MIP presolver workspace, which
793
-- the parameter ipp points to, freeing all the memory allocated to this
796
void ipp_delete_wksp(IPP *ipp)
797
{ if (ipp->row_pool != NULL) dmp_delete_pool(ipp->row_pool);
798
if (ipp->col_pool != NULL) dmp_delete_pool(ipp->col_pool);
799
if (ipp->aij_pool != NULL) dmp_delete_pool(ipp->aij_pool);
800
if (ipp->tqe_pool != NULL) dmp_delete_pool(ipp->tqe_pool);
801
if (ipp->col_stat != NULL) ufree(ipp->col_stat);
802
if (ipp->col_mipx != NULL) ufree(ipp->col_mipx);