diff options
Diffstat (limited to 'test/monniaux/glpk-4.65/src/npp/npp5.c')
| -rw-r--r-- | test/monniaux/glpk-4.65/src/npp/npp5.c | 809 |
1 files changed, 809 insertions, 0 deletions
diff --git a/test/monniaux/glpk-4.65/src/npp/npp5.c b/test/monniaux/glpk-4.65/src/npp/npp5.c new file mode 100644 index 00000000..2fad496d --- /dev/null +++ b/test/monniaux/glpk-4.65/src/npp/npp5.c @@ -0,0 +1,809 @@ +/* npp5.c */ + +/*********************************************************************** +* This code is part of GLPK (GNU Linear Programming Kit). +* +* Copyright (C) 2009-2017 Andrew Makhorin, Department for Applied +* Informatics, Moscow Aviation Institute, Moscow, Russia. All rights +* reserved. E-mail: <mao@gnu.org>. +* +* GLPK is free software: you can redistribute it and/or modify it +* under the terms of the GNU General Public License as published by +* the Free Software Foundation, either version 3 of the License, or +* (at your option) any later version. +* +* GLPK is distributed in the hope that it will be useful, but WITHOUT +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY +* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public +* License for more details. +* +* You should have received a copy of the GNU General Public License +* along with GLPK. If not, see <http://www.gnu.org/licenses/>. +***********************************************************************/ + +#include "env.h" +#include "npp.h" + +/*********************************************************************** +* NAME +* +* npp_clean_prob - perform initial LP/MIP processing +* +* SYNOPSIS +* +* #include "glpnpp.h" +* void npp_clean_prob(NPP *npp); +* +* DESCRIPTION +* +* The routine npp_clean_prob performs initial LP/MIP processing that +* currently includes: +* +* 1) removing free rows; +* +* 2) replacing double-sided constraint rows with almost identical +* bounds, by equality constraint rows; +* +* 3) removing fixed columns; +* +* 4) replacing double-bounded columns with almost identical bounds by +* fixed columns and removing those columns; +* +* 5) initial processing constraint coefficients (not implemented); +* +* 6) initial processing objective coefficients (not implemented). */ + +void npp_clean_prob(NPP *npp) +{ /* perform initial LP/MIP processing */ + NPPROW *row, *next_row; + NPPCOL *col, *next_col; + int ret; + xassert(npp == npp); + /* process rows which originally are free */ + for (row = npp->r_head; row != NULL; row = next_row) + { next_row = row->next; + if (row->lb == -DBL_MAX && row->ub == +DBL_MAX) + { /* process free row */ +#ifdef GLP_DEBUG + xprintf("1"); +#endif + npp_free_row(npp, row); + /* row was deleted */ + } + } + /* process rows which originally are double-sided inequalities */ + for (row = npp->r_head; row != NULL; row = next_row) + { next_row = row->next; + if (row->lb != -DBL_MAX && row->ub != +DBL_MAX && + row->lb < row->ub) + { ret = npp_make_equality(npp, row); + if (ret == 0) + ; + else if (ret == 1) + { /* row was replaced by equality constraint */ +#ifdef GLP_DEBUG + xprintf("2"); +#endif + } + else + xassert(ret != ret); + } + } + /* process columns which are originally fixed */ + for (col = npp->c_head; col != NULL; col = next_col) + { next_col = col->next; + if (col->lb == col->ub) + { /* process fixed column */ +#ifdef GLP_DEBUG + xprintf("3"); +#endif + npp_fixed_col(npp, col); + /* column was deleted */ + } + } + /* process columns which are originally double-bounded */ + for (col = npp->c_head; col != NULL; col = next_col) + { next_col = col->next; + if (col->lb != -DBL_MAX && col->ub != +DBL_MAX && + col->lb < col->ub) + { ret = npp_make_fixed(npp, col); + if (ret == 0) + ; + else if (ret == 1) + { /* column was replaced by fixed column; process it */ +#ifdef GLP_DEBUG + xprintf("4"); +#endif + npp_fixed_col(npp, col); + /* column was deleted */ + } + } + } + return; +} + +/*********************************************************************** +* NAME +* +* npp_process_row - perform basic row processing +* +* SYNOPSIS +* +* #include "glpnpp.h" +* int npp_process_row(NPP *npp, NPPROW *row, int hard); +* +* DESCRIPTION +* +* The routine npp_process_row performs basic row processing that +* currently includes: +* +* 1) removing empty row; +* +* 2) removing equality constraint row singleton and corresponding +* column; +* +* 3) removing inequality constraint row singleton and corresponding +* column if it was fixed; +* +* 4) performing general row analysis; +* +* 5) removing redundant row bounds; +* +* 6) removing forcing row and corresponding columns; +* +* 7) removing row which becomes free due to redundant bounds; +* +* 8) computing implied bounds for all columns in the row and using +* them to strengthen current column bounds (MIP only, optional, +* performed if the flag hard is on). +* +* Additionally the routine may activate affected rows and/or columns +* for further processing. +* +* RETURNS +* +* 0 success; +* +* GLP_ENOPFS primal/integer infeasibility detected; +* +* GLP_ENODFS dual infeasibility detected. */ + +int npp_process_row(NPP *npp, NPPROW *row, int hard) +{ /* perform basic row processing */ + NPPCOL *col; + NPPAIJ *aij, *next_aij, *aaa; + int ret; + /* row must not be free */ + xassert(!(row->lb == -DBL_MAX && row->ub == +DBL_MAX)); + /* start processing row */ + if (row->ptr == NULL) + { /* empty row */ + ret = npp_empty_row(npp, row); + if (ret == 0) + { /* row was deleted */ +#ifdef GLP_DEBUG + xprintf("A"); +#endif + return 0; + } + else if (ret == 1) + { /* primal infeasibility */ + return GLP_ENOPFS; + } + else + xassert(ret != ret); + } + if (row->ptr->r_next == NULL) + { /* row singleton */ + col = row->ptr->col; + if (row->lb == row->ub) + { /* equality constraint */ + ret = npp_eq_singlet(npp, row); + if (ret == 0) + { /* column was fixed, row was deleted */ +#ifdef GLP_DEBUG + xprintf("B"); +#endif + /* activate rows affected by column */ + for (aij = col->ptr; aij != NULL; aij = aij->c_next) + npp_activate_row(npp, aij->row); + /* process fixed column */ + npp_fixed_col(npp, col); + /* column was deleted */ + return 0; + } + else if (ret == 1 || ret == 2) + { /* primal/integer infeasibility */ + return GLP_ENOPFS; + } + else + xassert(ret != ret); + } + else + { /* inequality constraint */ + ret = npp_ineq_singlet(npp, row); + if (0 <= ret && ret <= 3) + { /* row was deleted */ +#ifdef GLP_DEBUG + xprintf("C"); +#endif + /* activate column, since its length was changed due to + row deletion */ + npp_activate_col(npp, col); + if (ret >= 2) + { /* column bounds changed significantly or column was + fixed */ + /* activate rows affected by column */ + for (aij = col->ptr; aij != NULL; aij = aij->c_next) + npp_activate_row(npp, aij->row); + } + if (ret == 3) + { /* column was fixed; process it */ +#ifdef GLP_DEBUG + xprintf("D"); +#endif + npp_fixed_col(npp, col); + /* column was deleted */ + } + return 0; + } + else if (ret == 4) + { /* primal infeasibility */ + return GLP_ENOPFS; + } + else + xassert(ret != ret); + } + } +#if 0 + /* sometimes this causes too large round-off errors; probably + pivot coefficient should be chosen more carefully */ + if (row->ptr->r_next->r_next == NULL) + { /* row doubleton */ + if (row->lb == row->ub) + { /* equality constraint */ + if (!(row->ptr->col->is_int || + row->ptr->r_next->col->is_int)) + { /* both columns are continuous */ + NPPCOL *q; + q = npp_eq_doublet(npp, row); + if (q != NULL) + { /* column q was eliminated */ +#ifdef GLP_DEBUG + xprintf("E"); +#endif + /* now column q is singleton of type "implied slack + variable"; we process it here to make sure that on + recovering basic solution the row is always active + equality constraint (as required by the routine + rcv_eq_doublet) */ + xassert(npp_process_col(npp, q) == 0); + /* column q was deleted; note that row p also may be + deleted */ + return 0; + } + } + } + } +#endif + /* general row analysis */ + ret = npp_analyze_row(npp, row); + xassert(0x00 <= ret && ret <= 0xFF); + if (ret == 0x33) + { /* row bounds are inconsistent with column bounds */ + return GLP_ENOPFS; + } + if ((ret & 0x0F) == 0x00) + { /* row lower bound does not exist or redundant */ + if (row->lb != -DBL_MAX) + { /* remove redundant row lower bound */ +#ifdef GLP_DEBUG + xprintf("F"); +#endif + npp_inactive_bound(npp, row, 0); + } + } + else if ((ret & 0x0F) == 0x01) + { /* row lower bound can be active */ + /* see below */ + } + else if ((ret & 0x0F) == 0x02) + { /* row lower bound is a forcing bound */ +#ifdef GLP_DEBUG + xprintf("G"); +#endif + /* process forcing row */ + if (npp_forcing_row(npp, row, 0) == 0) +fixup: { /* columns were fixed, row was made free */ + for (aij = row->ptr; aij != NULL; aij = next_aij) + { /* process column fixed by forcing row */ +#ifdef GLP_DEBUG + xprintf("H"); +#endif + col = aij->col; + next_aij = aij->r_next; + /* activate rows affected by column */ + for (aaa = col->ptr; aaa != NULL; aaa = aaa->c_next) + npp_activate_row(npp, aaa->row); + /* process fixed column */ + npp_fixed_col(npp, col); + /* column was deleted */ + } + /* process free row (which now is empty due to deletion of + all its columns) */ + npp_free_row(npp, row); + /* row was deleted */ + return 0; + } + } + else + xassert(ret != ret); + if ((ret & 0xF0) == 0x00) + { /* row upper bound does not exist or redundant */ + if (row->ub != +DBL_MAX) + { /* remove redundant row upper bound */ +#ifdef GLP_DEBUG + xprintf("I"); +#endif + npp_inactive_bound(npp, row, 1); + } + } + else if ((ret & 0xF0) == 0x10) + { /* row upper bound can be active */ + /* see below */ + } + else if ((ret & 0xF0) == 0x20) + { /* row upper bound is a forcing bound */ +#ifdef GLP_DEBUG + xprintf("J"); +#endif + /* process forcing row */ + if (npp_forcing_row(npp, row, 1) == 0) goto fixup; + } + else + xassert(ret != ret); + if (row->lb == -DBL_MAX && row->ub == +DBL_MAX) + { /* row became free due to redundant bounds removal */ +#ifdef GLP_DEBUG + xprintf("K"); +#endif + /* activate its columns, since their length will change due + to row deletion */ + for (aij = row->ptr; aij != NULL; aij = aij->r_next) + npp_activate_col(npp, aij->col); + /* process free row */ + npp_free_row(npp, row); + /* row was deleted */ + return 0; + } +#if 1 /* 23/XII-2009 */ + /* row lower and/or upper bounds can be active */ + if (npp->sol == GLP_MIP && hard) + { /* improve current column bounds (optional) */ + if (npp_improve_bounds(npp, row, 1) < 0) + return GLP_ENOPFS; + } +#endif + return 0; +} + +/*********************************************************************** +* NAME +* +* npp_improve_bounds - improve current column bounds +* +* SYNOPSIS +* +* #include "glpnpp.h" +* int npp_improve_bounds(NPP *npp, NPPROW *row, int flag); +* +* DESCRIPTION +* +* The routine npp_improve_bounds analyzes specified row (inequality +* or equality constraint) to determine implied column bounds and then +* uses these bounds to improve (strengthen) current column bounds. +* +* If the flag is on and current column bounds changed significantly +* or the column was fixed, the routine activate rows affected by the +* column for further processing. (This feature is intended to be used +* in the main loop of the routine npp_process_row.) +* +* NOTE: This operation can be used for MIP problem only. +* +* RETURNS +* +* The routine npp_improve_bounds returns the number of significantly +* changed bounds plus the number of column having been fixed due to +* bound improvements. However, if the routine detects primal/integer +* infeasibility, it returns a negative value. */ + +int npp_improve_bounds(NPP *npp, NPPROW *row, int flag) +{ /* improve current column bounds */ + NPPCOL *col; + NPPAIJ *aij, *next_aij, *aaa; + int kase, ret, count = 0; + double lb, ub; + xassert(npp->sol == GLP_MIP); + /* row must not be free */ + xassert(!(row->lb == -DBL_MAX && row->ub == +DBL_MAX)); + /* determine implied column bounds */ + npp_implied_bounds(npp, row); + /* and use these bounds to strengthen current column bounds */ + for (aij = row->ptr; aij != NULL; aij = next_aij) + { col = aij->col; + next_aij = aij->r_next; + for (kase = 0; kase <= 1; kase++) + { /* save current column bounds */ + lb = col->lb, ub = col->ub; + if (kase == 0) + { /* process implied column lower bound */ + if (col->ll.ll == -DBL_MAX) continue; + ret = npp_implied_lower(npp, col, col->ll.ll); + } + else + { /* process implied column upper bound */ + if (col->uu.uu == +DBL_MAX) continue; + ret = npp_implied_upper(npp, col, col->uu.uu); + } + if (ret == 0 || ret == 1) + { /* current column bounds did not change or changed, but + not significantly; restore current column bounds */ + col->lb = lb, col->ub = ub; + } + else if (ret == 2 || ret == 3) + { /* current column bounds changed significantly or column + was fixed */ +#ifdef GLP_DEBUG + xprintf("L"); +#endif + count++; + /* activate other rows affected by column, if required */ + if (flag) + { for (aaa = col->ptr; aaa != NULL; aaa = aaa->c_next) + { if (aaa->row != row) + npp_activate_row(npp, aaa->row); + } + } + if (ret == 3) + { /* process fixed column */ +#ifdef GLP_DEBUG + xprintf("M"); +#endif + npp_fixed_col(npp, col); + /* column was deleted */ + break; /* for kase */ + } + } + else if (ret == 4) + { /* primal/integer infeasibility */ + return -1; + } + else + xassert(ret != ret); + } + } + return count; +} + +/*********************************************************************** +* NAME +* +* npp_process_col - perform basic column processing +* +* SYNOPSIS +* +* #include "glpnpp.h" +* int npp_process_col(NPP *npp, NPPCOL *col); +* +* DESCRIPTION +* +* The routine npp_process_col performs basic column processing that +* currently includes: +* +* 1) fixing and removing empty column; +* +* 2) removing column singleton, which is implied slack variable, and +* corresponding row if it becomes free; +* +* 3) removing bounds of column, which is implied free variable, and +* replacing corresponding row by equality constraint. +* +* Additionally the routine may activate affected rows and/or columns +* for further processing. +* +* RETURNS +* +* 0 success; +* +* GLP_ENOPFS primal/integer infeasibility detected; +* +* GLP_ENODFS dual infeasibility detected. */ + +int npp_process_col(NPP *npp, NPPCOL *col) +{ /* perform basic column processing */ + NPPROW *row; + NPPAIJ *aij; + int ret; + /* column must not be fixed */ + xassert(col->lb < col->ub); + /* start processing column */ + if (col->ptr == NULL) + { /* empty column */ + ret = npp_empty_col(npp, col); + if (ret == 0) + { /* column was fixed and deleted */ +#ifdef GLP_DEBUG + xprintf("N"); +#endif + return 0; + } + else if (ret == 1) + { /* dual infeasibility */ + return GLP_ENODFS; + } + else + xassert(ret != ret); + } + if (col->ptr->c_next == NULL) + { /* column singleton */ + row = col->ptr->row; + if (row->lb == row->ub) + { /* equality constraint */ + if (!col->is_int) +slack: { /* implied slack variable */ +#ifdef GLP_DEBUG + xprintf("O"); +#endif + npp_implied_slack(npp, col); + /* column was deleted */ + if (row->lb == -DBL_MAX && row->ub == +DBL_MAX) + { /* row became free due to implied slack variable */ +#ifdef GLP_DEBUG + xprintf("P"); +#endif + /* activate columns affected by row */ + for (aij = row->ptr; aij != NULL; aij = aij->r_next) + npp_activate_col(npp, aij->col); + /* process free row */ + npp_free_row(npp, row); + /* row was deleted */ + } + else + { /* row became inequality constraint; activate it + since its length changed due to column deletion */ + npp_activate_row(npp, row); + } + return 0; + } + } + else + { /* inequality constraint */ + if (!col->is_int) + { ret = npp_implied_free(npp, col); + if (ret == 0) + { /* implied free variable */ +#ifdef GLP_DEBUG + xprintf("Q"); +#endif + /* column bounds were removed, row was replaced by + equality constraint */ + goto slack; + } + else if (ret == 1) + { /* column is not implied free variable, because its + lower and/or upper bounds can be active */ + } + else if (ret == 2) + { /* dual infeasibility */ + return GLP_ENODFS; + } + } + } + } + /* column still exists */ + return 0; +} + +/*********************************************************************** +* NAME +* +* npp_process_prob - perform basic LP/MIP processing +* +* SYNOPSIS +* +* #include "glpnpp.h" +* int npp_process_prob(NPP *npp, int hard); +* +* DESCRIPTION +* +* The routine npp_process_prob performs basic LP/MIP processing that +* currently includes: +* +* 1) initial LP/MIP processing (see the routine npp_clean_prob), +* +* 2) basic row processing (see the routine npp_process_row), and +* +* 3) basic column processing (see the routine npp_process_col). +* +* If the flag hard is on, the routine attempts to improve current +* column bounds multiple times within the main processing loop, in +* which case this feature may take a time. Otherwise, if the flag hard +* is off, improving column bounds is performed only once at the end of +* the main loop. (Note that this feature is used for MIP only.) +* +* The routine uses two sets: the set of active rows and the set of +* active columns. Rows/columns are marked by a flag (the field temp in +* NPPROW/NPPCOL). If the flag is non-zero, the row/column is active, +* in which case it is placed in the beginning of the row/column list; +* otherwise, if the flag is zero, the row/column is inactive, in which +* case it is placed in the end of the row/column list. If a row/column +* being currently processed may affect other rows/columns, the latters +* are activated for further processing. +* +* RETURNS +* +* 0 success; +* +* GLP_ENOPFS primal/integer infeasibility detected; +* +* GLP_ENODFS dual infeasibility detected. */ + +int npp_process_prob(NPP *npp, int hard) +{ /* perform basic LP/MIP processing */ + NPPROW *row; + NPPCOL *col; + int processing, ret; + /* perform initial LP/MIP processing */ + npp_clean_prob(npp); + /* activate all remaining rows and columns */ + for (row = npp->r_head; row != NULL; row = row->next) + row->temp = 1; + for (col = npp->c_head; col != NULL; col = col->next) + col->temp = 1; + /* main processing loop */ + processing = 1; + while (processing) + { processing = 0; + /* process all active rows */ + for (;;) + { row = npp->r_head; + if (row == NULL || !row->temp) break; + npp_deactivate_row(npp, row); + ret = npp_process_row(npp, row, hard); + if (ret != 0) goto done; + processing = 1; + } + /* process all active columns */ + for (;;) + { col = npp->c_head; + if (col == NULL || !col->temp) break; + npp_deactivate_col(npp, col); + ret = npp_process_col(npp, col); + if (ret != 0) goto done; + processing = 1; + } + } +#if 1 /* 23/XII-2009 */ + if (npp->sol == GLP_MIP && !hard) + { /* improve current column bounds (optional) */ + for (row = npp->r_head; row != NULL; row = row->next) + { if (npp_improve_bounds(npp, row, 0) < 0) + { ret = GLP_ENOPFS; + goto done; + } + } + } +#endif + /* all seems ok */ + ret = 0; +done: xassert(ret == 0 || ret == GLP_ENOPFS || ret == GLP_ENODFS); +#ifdef GLP_DEBUG + xprintf("\n"); +#endif + return ret; +} + +/**********************************************************************/ + +int npp_simplex(NPP *npp, const glp_smcp *parm) +{ /* process LP prior to applying primal/dual simplex method */ + int ret; + xassert(npp->sol == GLP_SOL); + xassert(parm == parm); + ret = npp_process_prob(npp, 0); + return ret; +} + +/**********************************************************************/ + +int npp_integer(NPP *npp, const glp_iocp *parm) +{ /* process MIP prior to applying branch-and-bound method */ + NPPROW *row, *prev_row; + NPPCOL *col; + NPPAIJ *aij; + int count, ret; + xassert(npp->sol == GLP_MIP); + xassert(parm == parm); + /*==============================================================*/ + /* perform basic MIP processing */ + ret = npp_process_prob(npp, 1); + if (ret != 0) goto done; + /*==============================================================*/ + /* binarize problem, if required */ + if (parm->binarize) + npp_binarize_prob(npp); + /*==============================================================*/ + /* identify hidden packing inequalities */ + count = 0; + /* new rows will be added to the end of the row list, so we go + from the end to beginning of the row list */ + for (row = npp->r_tail; row != NULL; row = prev_row) + { prev_row = row->prev; + /* skip free row */ + if (row->lb == -DBL_MAX && row->ub == +DBL_MAX) continue; + /* skip equality constraint */ + if (row->lb == row->ub) continue; + /* skip row having less than two variables */ + if (row->ptr == NULL || row->ptr->r_next == NULL) continue; + /* skip row having non-binary variables */ + for (aij = row->ptr; aij != NULL; aij = aij->r_next) + { col = aij->col; + if (!(col->is_int && col->lb == 0.0 && col->ub == 1.0)) + break; + } + if (aij != NULL) continue; + count += npp_hidden_packing(npp, row); + } + if (count > 0) + xprintf("%d hidden packing inequaliti(es) were detected\n", + count); + /*==============================================================*/ + /* identify hidden covering inequalities */ + count = 0; + /* new rows will be added to the end of the row list, so we go + from the end to beginning of the row list */ + for (row = npp->r_tail; row != NULL; row = prev_row) + { prev_row = row->prev; + /* skip free row */ + if (row->lb == -DBL_MAX && row->ub == +DBL_MAX) continue; + /* skip equality constraint */ + if (row->lb == row->ub) continue; + /* skip row having less than three variables */ + if (row->ptr == NULL || row->ptr->r_next == NULL || + row->ptr->r_next->r_next == NULL) continue; + /* skip row having non-binary variables */ + for (aij = row->ptr; aij != NULL; aij = aij->r_next) + { col = aij->col; + if (!(col->is_int && col->lb == 0.0 && col->ub == 1.0)) + break; + } + if (aij != NULL) continue; + count += npp_hidden_covering(npp, row); + } + if (count > 0) + xprintf("%d hidden covering inequaliti(es) were detected\n", + count); + /*==============================================================*/ + /* reduce inequality constraint coefficients */ + count = 0; + /* new rows will be added to the end of the row list, so we go + from the end to beginning of the row list */ + for (row = npp->r_tail; row != NULL; row = prev_row) + { prev_row = row->prev; + /* skip equality constraint */ + if (row->lb == row->ub) continue; + count += npp_reduce_ineq_coef(npp, row); + } + if (count > 0) + xprintf("%d constraint coefficient(s) were reduced\n", count); + /*==============================================================*/ +#ifdef GLP_DEBUG + routine(npp); +#endif + /*==============================================================*/ + /* all seems ok */ + ret = 0; +done: return ret; +} + +/* eof */ |