github.com/Interlisp.maiko
1
0
Fork 0
mirror of https://github.com/Interlisp/maiko.git synced 2026-01-14 15:36:34 +00:00
Code Issues Projects Releases Wiki Activity
Interlisp.maiko/src/lpkit.c
Nick Briggs 6fb4b9a189
Remove sccs id lines (#98) 
* Remove static char *id = from all source files.

The same information is included in a comment in each source file.

* Remove unused template file 'id'
2020-12-19 19:08:52 -08:00

1508 lines
44 KiB
C
Raw Blame History

/* $Id: lpkit.c,v 1.2 1999/01/03 02:07:18 sybalsky Exp $ (C) Copyright Venue, All Rights Reserved */
/************************************************************************/
/* */
/* (C) Copyright 1989-99 Venue. All Rights Reserved. */
/* Manufactured in the United States of America. */
/* */
/************************************************************************/
#include "version.h"
#include "lpkit.h"
#include "lpglob.h"
#include <stdlib.h>
#include <stdarg.h>
#include <stdio.h>
/* Globals */
lprec *Lp = NULL; /* pointer to active problem */
int Rows;
int Columns;
int Sum;
int Non_zeros;
int Level;
matrec *Mat;
int *Col_no;
int *Col_end;
int *Row_end;
REAL *Orig_rh;
REAL *Rh;
REAL *Rhs;
short *Must_be_int;
REAL *Orig_upbo;
REAL *Orig_lowbo;
REAL *Upbo;
REAL *Lowbo;
int *Bas;
short *Basis;
short *Lower;
int Eta_alloc;
int Eta_size;
REAL *Eta_value;
int *Eta_row_nr;
int *Eta_col_end;
int Num_inv;
REAL *Solution;
REAL *Best_solution;
REAL Infinite;
REAL Epsilon;
REAL Epsb;
REAL Epsd;
REAL Epsel;
REAL TREJ;
REAL TINV;
short Maximise;
short Floor_first;
REAL Extrad;
int Warn_count; /* used in CHECK version of rounding macro */
#ifdef NEVER
void error(char *format, ...) {
va_list ap;
va_start(ap, format);
vfprintf(stderr, format, ap);
fputc('\n', stderr);
va_end(ap);
exit(FAIL);
}
#endif
lprec *make_lp(int rows, int columns) {
lprec *newlp;
int i, sum;
sstate *st; /* to save state in MILP runs */
sum = rows + columns;
if (rows < 0 || columns < 0) error("rows < 0 or columns < 0");
CALLOC(newlp, 1, lprec);
strcpy(newlp->lp_name, "unnamed");
newlp->active = FALSE;
newlp->verbose = FALSE;
newlp->print_duals = FALSE;
newlp->print_sol = FALSE;
newlp->debug = FALSE;
newlp->print_at_invert = FALSE;
newlp->trace = FALSE;
newlp->rows = rows;
newlp->columns = columns;
newlp->sum = sum;
newlp->rows_alloc = rows;
newlp->columns_alloc = columns;
newlp->sum_alloc = sum;
newlp->names_used = FALSE;
newlp->obj_bound = DEF_INFINITE;
newlp->infinite = DEF_INFINITE;
newlp->epsilon = DEF_EPSILON;
newlp->epsb = DEF_EPSB;
newlp->epsd = DEF_EPSD;
newlp->epsel = DEF_EPSEL;
newlp->non_zeros = 0;
newlp->mat_alloc = 1;
CALLOC(newlp->mat, newlp->mat_alloc, matrec);
CALLOC(newlp->col_no, newlp->mat_alloc, int);
CALLOC(newlp->col_end, columns + 1, int);
CALLOC(newlp->row_end, rows + 1, int);
newlp->row_end_valid = FALSE;
CALLOC(newlp->orig_rh, rows + 1, REAL);
CALLOC(newlp->rh, rows + 1, REAL);
CALLOC(newlp->rhs, rows + 1, REAL);
CALLOC(newlp->must_be_int, sum + 1, short);
for (i = 0; i <= sum; i++) newlp->must_be_int[i] = FALSE;
CALLOC(newlp->orig_upbo, sum + 1, REAL);
for (i = 0; i <= sum; i++) newlp->orig_upbo[i] = newlp->infinite;
CALLOC(newlp->upbo, sum + 1, REAL);
CALLOC(newlp->orig_lowbo, sum + 1, REAL);
CALLOC(newlp->lowbo, sum + 1, REAL);
newlp->basis_valid = TRUE;
CALLOC(newlp->bas, rows + 1, int);
CALLOC(newlp->basis, sum + 1, short);
CALLOC(newlp->lower, sum + 1, short);
for (i = 0; i <= rows; i++) {
newlp->bas[i] = i;
newlp->basis[i] = TRUE;
}
for (i = rows + 1; i <= sum; i++) newlp->basis[i] = FALSE;
for (i = 0; i <= sum; i++) newlp->lower[i] = TRUE;
newlp->eta_valid = TRUE;
newlp->eta_size = 0;
newlp->eta_alloc = 10000;
newlp->max_num_inv = DEFNUMINV;
newlp->nr_lagrange = 0;
CALLOC(newlp->eta_value, newlp->eta_alloc, REAL);
CALLOC(newlp->eta_row_nr, newlp->eta_alloc, int);
CALLOC(newlp->eta_col_end, newlp->rows_alloc + newlp->max_num_inv, int);
newlp->bb_rule = FIRST_NI;
newlp->break_at_int = FALSE;
newlp->break_value = 0;
newlp->iter = 0;
newlp->total_iter = 0;
CALLOC(newlp->solution, sum + 1, REAL);
CALLOC(newlp->best_solution, sum + 1, REAL);
CALLOC(newlp->duals, rows + 1, REAL);
newlp->maximise = FALSE;
newlp->floor_first = TRUE;
newlp->scaling_used = FALSE;
newlp->columns_scaled = FALSE;
CALLOC(newlp->ch_sign, rows + 1, short);
for (i = 0; i <= rows; i++) newlp->ch_sign[i] = FALSE;
st = newlp->solve_states = NULL;
for (i = 0; i < Sum + 2; i++) {
if (st) /* There's an old state-saver to use; re-use it. */
{
st->saved = 0; /* clear the in-use field */
st->notint = 0; /* And the not-integer field */
st = st->next;
} else {
st = (sstate *)malloc(sizeof(sstate));
if (!st) ERROR(ERR_NOMEM); /* Tell the guy there's no memory */
st->next = newlp->solve_states;
st->saved = 0;
st->notint = 0; /* And the not-integer field */
newlp->solve_states = st;
st = (sstate *)NULL;
}
}
newlp->valid = FALSE;
return (newlp);
}
void delete_lp(lprec *lp) {
int i;
if (lp->active) Lp = NULL;
if (lp->names_used) {
free(lp->row_name);
free(lp->col_name);
}
free(lp->mat);
free(lp->col_no);
free(lp->col_end);
free(lp->row_end);
free(lp->orig_rh);
free(lp->rh);
free(lp->rhs);
free(lp->must_be_int);
free(lp->orig_upbo);
free(lp->orig_lowbo);
free(lp->upbo);
free(lp->lowbo);
free(lp->bas);
free(lp->basis);
free(lp->lower);
free(lp->eta_value);
free(lp->eta_row_nr);
free(lp->eta_col_end);
free(lp->solution);
free(lp->best_solution);
free(lp->duals);
free(lp->ch_sign);
if (lp->scaling_used) free(lp->scale);
if (lp->nr_lagrange > 0) {
free(lp->lag_rhs);
free(lp->lambda);
free(lp->lag_con_type);
for (i = 0; i < lp->nr_lagrange; i++) free(lp->lag_row[i]);
free(lp->lag_row);
}
free(lp);
lp = NULL;
}
lprec *copy_lp(lprec *lp) {
lprec *newlp;
int i, rowsplus, colsplus, sumplus;
rowsplus = lp->rows_alloc + 1;
colsplus = lp->columns_alloc + 1;
sumplus = lp->sum_alloc + 1;
MALLOCCPY(newlp, lp, 1, lprec); /* copy all non pointers */
newlp->active = FALSE;
if (newlp->names_used) {
MALLOCCPY(newlp->col_name, lp->col_name, colsplus, nstring);
MALLOCCPY(newlp->row_name, lp->row_name, rowsplus, nstring);
}
MALLOCCPY(newlp->mat, lp->mat, newlp->mat_alloc, matrec);
MALLOCCPY(newlp->col_end, lp->col_end, colsplus, int);
MALLOCCPY(newlp->col_no, lp->col_no, newlp->mat_alloc, int);
MALLOCCPY(newlp->row_end, lp->row_end, rowsplus, int);
MALLOCCPY(newlp->orig_rh, lp->orig_rh, rowsplus, REAL);
MALLOCCPY(newlp->rh, lp->rh, rowsplus, REAL);
MALLOCCPY(newlp->rhs, lp->rhs, rowsplus, REAL);
MALLOCCPY(newlp->must_be_int, lp->must_be_int, sumplus, short);
MALLOCCPY(newlp->orig_upbo, lp->orig_upbo, sumplus, REAL);
MALLOCCPY(newlp->orig_lowbo, lp->orig_lowbo, sumplus, REAL);
MALLOCCPY(newlp->upbo, lp->upbo, sumplus, REAL);
MALLOCCPY(newlp->lowbo, lp->lowbo, sumplus, REAL);
MALLOCCPY(newlp->bas, lp->bas, rowsplus, int);
MALLOCCPY(newlp->basis, lp->basis, sumplus, short);
MALLOCCPY(newlp->lower, lp->lower, sumplus, short);
MALLOCCPY(newlp->eta_value, lp->eta_value, lp->eta_alloc, REAL);
MALLOCCPY(newlp->eta_row_nr, lp->eta_row_nr, lp->eta_alloc, int);
MALLOCCPY(newlp->eta_col_end, lp->eta_col_end, lp->rows_alloc + lp->max_num_inv, int);
MALLOCCPY(newlp->solution, lp->solution, sumplus, REAL);
MALLOCCPY(newlp->best_solution, lp->best_solution, sumplus, REAL);
MALLOCCPY(newlp->duals, lp->duals, rowsplus, REAL);
MALLOCCPY(newlp->ch_sign, lp->ch_sign, rowsplus, short);
if (newlp->scaling_used) MALLOCCPY(newlp->scale, lp->scale, sumplus, REAL);
if (newlp->nr_lagrange > 0) {
MALLOCCPY(newlp->lag_rhs, lp->lag_rhs, newlp->nr_lagrange, REAL);
MALLOCCPY(newlp->lambda, lp->lambda, newlp->nr_lagrange, REAL);
MALLOCCPY(newlp->lag_con_type, lp->lag_con_type, newlp->nr_lagrange, short);
MALLOC(newlp->lag_row, newlp->nr_lagrange, REAL *);
for (i = 0; i < newlp->nr_lagrange; i++)
MALLOCCPY(newlp->lag_row[i], lp->lag_row[i], colsplus, REAL);
}
return (newlp);
}
void inc_mat_space(lprec *lp, int maxextra) {
if (lp->non_zeros + maxextra >= lp->mat_alloc) {
lp->mat_alloc = lp->non_zeros + maxextra;
REALLOC(lp->mat, lp->mat_alloc, matrec);
REALLOC(lp->col_no, lp->mat_alloc, int);
if (lp->active) {
Mat = lp->mat;
Col_no = lp->col_no;
}
}
}
void inc_row_space(lprec *lp) {
if (lp->rows > lp->rows_alloc) {
lp->rows_alloc = lp->rows + 10;
lp->sum_alloc = lp->rows_alloc + lp->columns_alloc;
REALLOC(lp->orig_rh, lp->rows_alloc + 1, REAL);
REALLOC(lp->rh, lp->rows_alloc + 1, REAL);
REALLOC(lp->rhs, lp->rows_alloc + 1, REAL);
REALLOC(lp->orig_upbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->upbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->orig_lowbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->lowbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->solution, lp->sum_alloc + 1, REAL);
REALLOC(lp->best_solution, lp->sum_alloc + 1, REAL);
REALLOC(lp->row_end, lp->rows_alloc + 1, int);
REALLOC(lp->basis, lp->sum_alloc + 1, short);
REALLOC(lp->lower, lp->sum_alloc + 1, short);
REALLOC(lp->must_be_int, lp->sum_alloc + 1, short);
REALLOC(lp->bas, lp->rows_alloc + 1, int);
REALLOC(lp->duals, lp->rows_alloc + 1, REAL);
REALLOC(lp->ch_sign, lp->rows_alloc + 1, short);
REALLOC(lp->eta_col_end, lp->rows_alloc + lp->max_num_inv, int);
if (lp->names_used) REALLOC(lp->row_name, lp->rows_alloc + 1, nstring);
if (lp->scaling_used) REALLOC(lp->scale, lp->sum_alloc + 1, REAL);
if (lp->active) set_globals(lp);
}
}
void inc_col_space(lprec *lp) {
if (lp->columns >= lp->columns_alloc) {
lp->columns_alloc = lp->columns + 10;
lp->sum_alloc = lp->rows_alloc + lp->columns_alloc;
REALLOC(lp->must_be_int, lp->sum_alloc + 1, short);
REALLOC(lp->orig_upbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->upbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->orig_lowbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->lowbo, lp->sum_alloc + 1, REAL);
REALLOC(lp->solution, lp->sum_alloc + 1, REAL);
REALLOC(lp->best_solution, lp->sum_alloc + 1, REAL);
REALLOC(lp->basis, lp->sum_alloc + 1, short);
REALLOC(lp->lower, lp->sum_alloc + 1, short);
if (lp->names_used) REALLOC(lp->col_name, lp->columns_alloc + 1, nstring);
if (lp->scaling_used) REALLOC(lp->scale, lp->sum_alloc + 1, REAL);
REALLOC(lp->col_end, lp->columns_alloc + 1, int);
if (lp->active) set_globals(lp);
}
}
void set_mat(lprec *lp, int Row, int Column, REAL Value) {
int elmnr, lastelm, i;
if (Row > lp->rows || Row < 0) error("Row out of range");
if (Column > lp->columns || Column < 1) error("Column out of range");
if (lp->scaling_used) Value *= lp->scale[Row] * lp->scale[lp->rows + Column];
if (TRUE /*abs(Value) > lp->epsilon*/) {
if (lp->basis[Column] == TRUE && Row > 0) lp->basis_valid = FALSE;
lp->eta_valid = FALSE;
elmnr = lp->col_end[Column - 1];
while ((elmnr < lp->col_end[Column]) ? (lp->mat[elmnr].row_nr != Row) : FALSE) elmnr++;
if ((elmnr != lp->col_end[Column]) ? (lp->mat[elmnr].row_nr == Row) : FALSE)
if (lp->scaling_used) {
if (lp->ch_sign[Row])
lp->mat[elmnr].value = -Value * lp->scale[Row] * lp->scale[Column];
else
lp->mat[elmnr].value = Value * lp->scale[Row] * lp->scale[Column];
} else {
if (lp->ch_sign[Row])
lp->mat[elmnr].value = -Value;
else
lp->mat[elmnr].value = Value;
}
else {
/* check if more space is needed for matrix */
inc_mat_space(lp, 1);
/* Shift the matrix */
lastelm = lp->non_zeros;
for (i = lastelm; i > elmnr; i--) lp->mat[i] = lp->mat[i - 1];
for (i = Column; i <= lp->columns; i++) lp->col_end[i]++;
/* Set new element */
lp->mat[elmnr].row_nr = Row;
if (lp->scaling_used) {
if (lp->ch_sign[Row])
lp->mat[elmnr].value = -Value * lp->scale[Row] * lp->scale[Column];
else
lp->mat[elmnr].value = Value * lp->scale[Row] * lp->scale[Column];
} else {
if (lp->ch_sign[Row])
lp->mat[elmnr].value = -Value;
else
lp->mat[elmnr].value = Value;
}
lp->row_end_valid = FALSE;
lp->non_zeros++;
if (lp->active) Non_zeros = lp->non_zeros;
}
}
}
void set_obj_fn(lprec *lp, REAL *row) {
int i;
for (i = 1; i <= lp->columns; i++) set_mat(lp, 0, i, row[i]);
}
void str_set_obj_fn(lprec *lp, char *row) {
int i;
REAL *arow;
char *p, *newp;
CALLOC(arow, lp->columns + 1, REAL);
p = row;
for (i = 1; i <= lp->columns; i++) {
arow[i] = (REAL)strtod(p, &newp);
if (p == newp)
error("Bad string in str_set_obj_fn");
else
p = newp;
}
set_obj_fn(lp, arow);
free(arow);
}
void add_constraint(lprec *lp, REAL *row, short constr_type, REAL rh) {
matrec *newmat;
int i, j;
int elmnr;
int stcol;
int *addtoo;
int *num, *rownum, row_nr;
REAL theta;
MALLOC(addtoo, lp->columns + 1, int)
for (i = 1; i <= lp->columns; i++)
if (row[i] != 0) {
addtoo[i] = TRUE;
lp->non_zeros++;
} else
addtoo[i] = FALSE;
MALLOC(newmat, lp->non_zeros, matrec);
inc_mat_space(lp, 0);
lp->rows++;
lp->sum++;
inc_row_space(lp);
if (lp->scaling_used) {
/* shift scale */
for (i = lp->sum; i > lp->rows; i--) lp->scale[i] = lp->scale[i - 1];
lp->scale[lp->rows] = 1;
}
if (lp->names_used) sprintf(lp->row_name[lp->rows], "r_%d", lp->rows);
if (lp->scaling_used && lp->columns_scaled)
for (i = 1; i <= lp->columns; i++) row[i] *= lp->scale[lp->rows + i];
if (constr_type == GE)
lp->ch_sign[lp->rows] = TRUE;
else
lp->ch_sign[lp->rows] = FALSE;
elmnr = 0;
stcol = 0;
for (i = 1; i <= lp->columns; i++) {
for (j = stcol; j < lp->col_end[i]; j++) {
newmat[elmnr].row_nr = lp->mat[j].row_nr;
newmat[elmnr].value = lp->mat[j].value;
elmnr++;
}
if (addtoo[i]) {
if (lp->ch_sign[lp->rows])
newmat[elmnr].value = -row[i];
else
newmat[elmnr].value = row[i];
newmat[elmnr].row_nr = lp->rows;
elmnr++;
}
stcol = lp->col_end[i];
lp->col_end[i] = elmnr;
}
memcpy(lp->mat, newmat, lp->non_zeros * sizeof(matrec));
free(newmat);
free(addtoo);
for (i = lp->sum; i > lp->rows; i--) {
lp->orig_upbo[i] = lp->orig_upbo[i - 1];
lp->orig_lowbo[i] = lp->orig_lowbo[i - 1];
lp->basis[i] = lp->basis[i - 1];
lp->lower[i] = lp->lower[i - 1];
lp->must_be_int[i] = lp->must_be_int[i - 1];
}
for (i = 1; i <= lp->rows; i++)
if (lp->bas[i] >= lp->rows) lp->bas[i]++;
if (constr_type == LE || constr_type == GE) {
lp->orig_upbo[lp->rows] = lp->infinite;
} else if (constr_type == EQ) {
lp->orig_upbo[lp->rows] = 0;
} else {
fprintf(stderr, "Wrong constraint type\n");
exit(FAIL);
}
lp->orig_lowbo[lp->rows] = 0;
if (constr_type == GE && rh != 0)
lp->orig_rh[lp->rows] = -rh;
else
lp->orig_rh[lp->rows] = rh;
lp->row_end_valid = FALSE;
lp->bas[lp->rows] = lp->rows;
lp->basis[lp->rows] = TRUE;
lp->lower[lp->rows] = TRUE;
if (lp->active) set_globals(lp);
lp->eta_valid = FALSE;
}
void str_add_constraint(lprec *lp, char *row_string, short constr_type, REAL rh) {
int i;
REAL *aRow;
char *p, *newp;
CALLOC(aRow, lp->columns + 1, REAL);
p = row_string;
for (i = 1; i <= lp->columns; i++) {
aRow[i] = (REAL)strtod(p, &newp);
if (p == newp)
error("Bad string in str_add_constr");
else
p = newp;
}
add_constraint(lp, aRow, constr_type, rh);
free(aRow);
}
void del_constraint(lprec *lp, int del_row) {
int i, j;
unsigned elmnr;
int startcol;
int row_nr;
if (del_row < 1 || del_row > lp->rows) {
fprintf(stderr, "There is no constraint nr. %d\n", del_row);
exit(FAIL);
}
elmnr = 0;
startcol = 0;
for (i = 1; i <= lp->columns; i++) {
for (j = startcol; j < lp->col_end[i]; j++) {
if (lp->mat[j].row_nr != del_row) {
lp->mat[elmnr] = lp->mat[j];
if (lp->mat[elmnr].row_nr > del_row) lp->mat[elmnr].row_nr--;
elmnr++;
} else
lp->non_zeros--;
}
startcol = lp->col_end[i];
lp->col_end[i] = elmnr;
}
for (i = del_row; i < lp->rows; i++) {
lp->orig_rh[i] = lp->orig_rh[i + 1];
lp->ch_sign[i] = lp->ch_sign[i + 1];
lp->bas[i] = lp->bas[i + 1];
if (lp->names_used) strcpy(lp->row_name[i], lp->row_name[i + 1]);
}
for (i = 1; i < lp->rows; i++)
if (lp->bas[i] > del_row) lp->bas[i]--;
for (i = del_row; i < lp->sum; i++) {
lp->lower[i] = lp->lower[i + 1];
lp->basis[i] = lp->basis[i + 1];
lp->orig_upbo[i] = lp->orig_upbo[i + 1];
lp->orig_lowbo[i] = lp->orig_lowbo[i + 1];
lp->must_be_int[i] = lp->must_be_int[i + 1];
if (lp->scaling_used) lp->scale[i] = lp->scale[i + 1];
}
lp->rows--;
lp->sum--;
lp->row_end_valid = FALSE;
if (lp->active) set_globals(lp);
lp->eta_valid = FALSE;
lp->basis_valid = FALSE;
}
void add_lag_con(lprec *lp, REAL *row, short con_type, REAL rhs) {
int i;
REAL sign;
if (con_type == LE || con_type == EQ)
sign = 1;
else if (con_type == GE)
sign = -1;
else
error("con_type not implemented\n");
lp->nr_lagrange++;
if (lp->nr_lagrange == 1) {
CALLOC(lp->lag_row, lp->nr_lagrange, REAL *);
CALLOC(lp->lag_rhs, lp->nr_lagrange, REAL);
CALLOC(lp->lambda, lp->nr_lagrange, REAL);
CALLOC(lp->lag_con_type, lp->nr_lagrange, short);
} else {
REALLOC(lp->lag_row, lp->nr_lagrange, REAL *);
REALLOC(lp->lag_rhs, lp->nr_lagrange, REAL);
REALLOC(lp->lambda, lp->nr_lagrange, REAL);
REALLOC(lp->lag_con_type, lp->nr_lagrange, short);
}
CALLOC(lp->lag_row[lp->nr_lagrange - 1], lp->columns + 1, REAL);
lp->lag_rhs[lp->nr_lagrange - 1] = rhs * sign;
for (i = 1; i <= lp->columns; i++) lp->lag_row[lp->nr_lagrange - 1][i] = row[i] * sign;
lp->lambda[lp->nr_lagrange - 1] = 0;
lp->lag_con_type[lp->nr_lagrange - 1] = (con_type == EQ);
}
void str_add_lag_con(lprec *lp, char *row, short con_type, REAL rhs) {
int i;
REAL *a_row;
char *p, *new_p;
CALLOC(a_row, lp->columns + 1, REAL);
p = row;
for (i = 1; i <= lp->columns; i++) {
a_row[i] = (REAL)strtod(p, &new_p);
if (p == new_p)
error("Bad string in str_add_lag_con");
else
p = new_p;
}
add_lag_con(lp, a_row, con_type, rhs);
free(a_row);
}
void add_column(lprec *lp, REAL *column) {
int i, j, elmnr, row_nr, *num, *rownum;
lp->columns++;
lp->sum++;
inc_col_space(lp);
inc_mat_space(lp, lp->rows + 1);
if (lp->scaling_used) {
for (i = 0; i <= lp->rows; i++) column[i] *= lp->scale[i];
lp->scale[lp->sum] = 1;
}
elmnr = lp->col_end[lp->columns - 1];
for (i = 0; i <= lp->rows; i++)
if (column[i] != 0) {
lp->mat[elmnr].row_nr = i;
if (lp->ch_sign[i])
lp->mat[elmnr].value = -column[i];
else
lp->mat[elmnr].value = column[i];
lp->non_zeros++;
elmnr++;
}
lp->col_end[lp->columns] = elmnr;
lp->orig_lowbo[lp->sum] = 0;
lp->orig_upbo[lp->sum] = lp->infinite;
lp->lower[lp->sum] = TRUE;
lp->basis[lp->sum] = FALSE;
lp->must_be_int[lp->sum] = FALSE;
if (lp->names_used) sprintf(lp->col_name[lp->columns], "var_%d", lp->columns);
lp->row_end_valid = FALSE;
if (lp->active) {
Sum = lp->sum;
Columns = lp->columns;
Non_zeros = lp->non_zeros;
}
}
void str_add_column(lprec *lp, char *col_string) {
int i;
REAL *aCol;
char *p, *newp;
CALLOC(aCol, lp->rows + 1, REAL);
p = col_string;
for (i = 0; i <= lp->rows; i++) {
aCol[i] = (REAL)strtod(p, &newp);
if (p == newp)
error("Bad string in str_add_column");
else
p = newp;
}
add_column(lp, aCol);
free(aCol);
}
void del_column(lprec *lp, int column) {
int i, j, from_elm, to_elm, elm_in_col;
if (column > lp->columns || column < 1) error("Column out of range in del_column");
for (i = 1; i <= lp->rows; i++) {
if (lp->bas[i] == lp->rows + column)
lp->basis_valid = FALSE;
else if (lp->bas[i] > lp->rows + column)
lp->bas[i]--;
}
for (i = lp->rows + column; i < lp->sum; i++) {
if (lp->names_used) strcpy(lp->col_name[i - lp->rows], lp->col_name[i - lp->rows + 1]);
lp->must_be_int[i] = lp->must_be_int[i + 1];
lp->orig_upbo[i] = lp->orig_upbo[i + 1];
lp->orig_lowbo[i] = lp->orig_lowbo[i + 1];
lp->upbo[i] = lp->upbo[i + 1];
lp->lowbo[i] = lp->lowbo[i + 1];
lp->basis[i] = lp->basis[i + 1];
lp->lower[i] = lp->lower[i + 1];
if (lp->scaling_used) lp->scale[i] = lp->scale[i + 1];
}
for (i = 0; i < lp->nr_lagrange; i++)
for (j = column; j <= lp->columns; j++) lp->lag_row[i][j] = lp->lag_row[i][j + 1];
to_elm = lp->col_end[column - 1];
from_elm = lp->col_end[column];
elm_in_col = from_elm - to_elm;
for (i = from_elm; i < lp->non_zeros; i++) {
lp->mat[to_elm] = lp->mat[i];
to_elm++;
}
for (i = column; i < lp->columns; i++) lp->col_end[i] = lp->col_end[i + 1] - elm_in_col;
lp->non_zeros -= elm_in_col;
lp->row_end_valid = FALSE;
lp->eta_valid = FALSE;
lp->sum--;
lp->columns--;
if (lp->active) set_globals(lp);
}
void set_upbo(lprec *lp, int column, REAL value) {
if (column > lp->columns || column < 1) error("Column out of range");
if (lp->scaling_used) value /= lp->scale[lp->rows + column];
if (value < lp->orig_lowbo[lp->rows + column]) error("Upperbound must be >= lowerbound");
lp->eta_valid = FALSE;
lp->orig_upbo[lp->rows + column] = value;
}
void set_lowbo(lprec *lp, int column, REAL value) {
if (column > lp->columns || column < 1) error("Column out of range");
if (lp->scaling_used) value /= lp->scale[lp->rows + column];
if (value > lp->orig_upbo[lp->rows + column]) error("Upperbound must be >= lowerbound");
lp->eta_valid = FALSE;
lp->orig_lowbo[lp->rows + column] = value;
}
void set_int(lprec *lp, int column, short must_be_int) {
if (column > lp->columns || column < 1) error("Column out of range");
lp->must_be_int[lp->rows + column] = must_be_int;
if (must_be_int == TRUE)
if (lp->columns_scaled) unscale_columns(lp);
}
void set_rh(lprec *lp, int row, REAL value) {
if (row > lp->rows || row < 0) error("Row out of Range");
if (row == 0) /* setting of RHS of OF not meaningful */
{
fprintf(stderr, "Warning: attempt to set RHS of objective function, ignored\n");
return;
}
if (lp->scaling_used)
if (lp->ch_sign[row])
lp->orig_rh[row] = -value * lp->scale[row];
else
lp->orig_rh[row] = value * lp->scale[row];
else if (lp->ch_sign[row])
lp->orig_rh[row] = -value;
else
lp->orig_rh[row] = value;
lp->eta_valid = FALSE;
}
void set_rh_vec(lprec *lp, REAL *rh) {
int i;
if (lp->scaling_used)
for (i = 1; i <= lp->rows; i++)
if (lp->ch_sign[i])
lp->orig_rh[i] = -rh[i] * lp->scale[i];
else
lp->orig_rh[i] = rh[i] * lp->scale[i];
else
for (i = 1; i <= lp->rows; i++)
if (lp->ch_sign[i])
lp->orig_rh[i] = -rh[i];
else
lp->orig_rh[i] = rh[i];
lp->eta_valid = FALSE;
}
void str_set_rh_vec(lprec *lp, char *rh_string) {
int i;
REAL *newrh;
char *p, *newp;
CALLOC(newrh, lp->rows + 1, REAL);
p = rh_string;
for (i = 1; i <= lp->rows; i++) {
newrh[i] = (REAL)strtod(p, &newp);
if (p == newp)
error("Bad string in str_set_rh_vec");
else
p = newp;
}
set_rh_vec(lp, newrh);
free(newrh);
}
void set_maxim(lprec *lp) {
int i;
if (lp->maximise == FALSE) {
for (i = 0; i < lp->non_zeros; i++)
if (lp->mat[i].row_nr == 0) lp->mat[i].value *= -1;
lp->eta_valid = FALSE;
}
lp->maximise = TRUE;
lp->ch_sign[0] = TRUE;
if (lp->active) Maximise = TRUE;
}
void set_minim(lprec *lp) {
int i;
if (lp->maximise == TRUE) {
for (i = 0; i < lp->non_zeros; i++)
if (lp->mat[i].row_nr == 0) lp->mat[i].value = -lp->mat[i].value;
lp->eta_valid = FALSE;
}
lp->maximise = FALSE;
lp->ch_sign[0] = FALSE;
if (lp->active) Maximise = FALSE;
}
void set_constr_type(lprec *lp, int row, short con_type) {
int i;
if (row > lp->rows || row < 1) error("Row out of Range");
if (con_type == EQ) {
lp->orig_upbo[row] = 0;
lp->basis_valid = FALSE;
if (lp->ch_sign[row]) {
for (i = 0; i < lp->non_zeros; i++)
if (lp->mat[i].row_nr == row) lp->mat[i].value *= -1;
lp->eta_valid = FALSE;
lp->ch_sign[row] = FALSE;
if (lp->orig_rh[row] != 0) lp->orig_rh[row] *= -1;
}
} else if (con_type == LE) {
lp->orig_upbo[row] = lp->infinite;
lp->basis_valid = FALSE;
if (lp->ch_sign[row]) {
for (i = 0; i < lp->non_zeros; i++)
if (lp->mat[i].row_nr == row) lp->mat[i].value *= -1;
lp->eta_valid = FALSE;
lp->ch_sign[row] = FALSE;
if (lp->orig_rh[row] != 0) lp->orig_rh[row] *= -1;
}
} else if (con_type == GE) {
lp->orig_upbo[row] = lp->infinite;
lp->basis_valid = FALSE;
if (!lp->ch_sign[row]) {
for (i = 0; i < lp->non_zeros; i++)
if (lp->mat[i].row_nr == row) lp->mat[i].value *= -1;
lp->eta_valid = FALSE;
lp->ch_sign[row] = TRUE;
if (lp->orig_rh[row] != 0) lp->orig_rh[row] *= -1;
}
} else
error("Constraint type not (yet) implemented");
}
REAL mat_elm(lprec *lp, int row, int column) {
REAL value;
int elmnr;
if (row < 0 || row > lp->rows) error("Row out of range in mat_elm");
if (column < 1 || column > lp->columns) error("Column out of range in mat_elm");
value = 0;
elmnr = lp->col_end[column - 1];
while (lp->mat[elmnr].row_nr != row && elmnr < lp->col_end[column]) elmnr++;
if (elmnr != lp->col_end[column]) {
value = lp->mat[elmnr].value;
if (lp->ch_sign[row]) value = -value;
if (lp->scaling_used) value /= lp->scale[row] * lp->scale[lp->rows + column];
}
return (value);
}
void get_row(lprec *lp, int row_nr, REAL *row) {
int i, j;
REAL RowMult;
if (row_nr < 0 || row_nr > lp->rows) error("Row nr. out of range in get_row");
for (i = 1; i <= lp->columns; i++) {
row[i] = 0;
for (j = lp->col_end[i - 1]; j < lp->col_end[i]; j++)
if (lp->mat[j].row_nr == row_nr) row[i] = lp->mat[j].value;
if (lp->scaling_used) row[i] /= lp->scale[lp->rows + i] * lp->scale[row_nr];
}
if (lp->ch_sign[row_nr])
for (i = 0; i <= lp->columns; i++)
if (row[i] != 0) row[i] = -row[i];
}
void get_column(lprec *lp, int col_nr, REAL *column) {
int i;
if (col_nr < 1 || col_nr > lp->columns) error("Col. nr. out of range in get_column");
for (i = 0; i <= lp->rows; i++) column[i] = 0;
for (i = lp->col_end[col_nr - 1]; i < lp->col_end[col_nr]; i++)
column[lp->mat[i].row_nr] = lp->mat[i].value;
for (i = 0; i <= lp->rows; i++)
if (column[i] != 0) {
if (lp->ch_sign[i]) column[i] *= -1;
if (lp->scaling_used) column[i] /= (lp->scale[i] * lp->scale[lp->rows + col_nr]);
}
}
void get_reduced_costs(lprec *lp, REAL *rc) {
int varnr, i, j;
REAL f;
if (!lp->basis_valid) error("Not a valid basis in get_reduced_costs");
set_globals(lp);
if (!lp->eta_valid) invert();
for (i = 1; i <= lp->sum; i++) rc[i] = 0;
rc[0] = 1;
btran(rc);
for (i = 1; i <= lp->columns; i++) {
varnr = lp->rows + i;
if (!Basis[varnr])
if (Upbo[varnr] > 0) {
f = 0;
for (j = Col_end[i - 1]; j < Col_end[i]; j++) f += rc[Mat[j].row_nr] * Mat[j].value;
rc[varnr] = f;
}
}
for (i = 1; i <= Sum; i++) my_round(rc[i], Epsd);
}
short is_feasible(lprec *lp, REAL *values) {
int i, elmnr;
REAL *this_rhs;
REAL dist;
if (lp->scaling_used) {
for (i = lp->rows + 1; i <= lp->sum; i++)
if (values[i - lp->rows] < lp->orig_lowbo[i] * lp->scale[i] ||
values[i - lp->rows] > lp->orig_upbo[i] * lp->scale[i])
return (FALSE);
} else {
for (i = lp->rows + 1; i <= lp->sum; i++)
if (values[i - lp->rows] < lp->orig_lowbo[i] || values[i - lp->rows] > lp->orig_upbo[i])
return (FALSE);
}
CALLOC(this_rhs, lp->rows + 1, REAL)
for (i = 1; i <= lp->columns; i++)
for (elmnr = lp->col_end[i - 1]; elmnr < lp->col_end[i]; elmnr++)
this_rhs[lp->mat[elmnr].row_nr] += lp->mat[elmnr].value * values[i];
for (i = 1; i <= lp->rows; i++) {
dist = lp->orig_rh[i] - this_rhs[i];
my_round(dist, 0.001) if ((lp->orig_upbo[i] == 0 && dist != 0) || dist < 0) {
free(this_rhs);
return (FALSE);
}
}
free(this_rhs);
return (TRUE);
}
short column_in_lp(lprec *lp, REAL *testcolumn) {
int i, j;
short ident;
REAL value;
if (lp->scaling_used)
for (i = 1; i <= lp->columns; i++) {
ident = TRUE;
j = lp->col_end[i - 1];
while (ident && (j < lp->col_end[i])) {
value = lp->mat[j].value;
if (lp->ch_sign[lp->mat[j].row_nr]) value = -value;
value /= lp->scale[lp->rows + i];
value /= lp->scale[lp->mat[j].row_nr];
value -= testcolumn[lp->mat[j].row_nr];
if (my_abs(value) > 0.001) /* should be some epsilon? MB */
ident = FALSE;
j++;
}
if (ident) return (TRUE);
}
else
for (i = 1; i <= lp->columns; i++) {
ident = TRUE;
j = lp->col_end[i - 1];
while (ident && j < lp->col_end[i]) {
value = lp->mat[j].value;
if (lp->ch_sign[lp->mat[j].row_nr]) value *= -1;
value -= testcolumn[lp->mat[j].row_nr];
if (my_abs(value) > 0.001) ident = FALSE;
j++;
}
if (ident) return (TRUE);
}
return (FALSE);
}
void print_lp(lprec *lp) {
int i, j;
REAL *fatmat;
CALLOC(fatmat, (lp->rows + 1) * lp->columns, REAL);
for (i = 1; i <= lp->columns; i++)
for (j = lp->col_end[i - 1]; j < lp->col_end[i]; j++)
fatmat[(i - 1) * (lp->rows + 1) + lp->mat[j].row_nr] = lp->mat[j].value;
printf("problem name: %s\n", lp->lp_name);
printf(" ");
for (j = 1; j <= lp->columns; j++)
if (lp->names_used)
printf("%8s ", lp->col_name[j]);
else
printf("Var[%3d] ", j);
if (lp->maximise) {
printf("\nMaximise ");
for (j = 0; j < lp->columns; j++) printf("% 8.2f ", -fatmat[j * (lp->rows + 1)]);
} else {
printf("\nMinimize ");
for (j = 0; j < lp->columns; j++) printf("% 8.2f ", fatmat[j * (lp->rows + 1)]);
}
printf("\n");
for (i = 1; i <= lp->rows; i++) {
if (lp->names_used)
printf("%9s ", lp->row_name[i]);
else
printf("Row[%3d] ", i);
for (j = 0; j < lp->columns; j++)
if (lp->ch_sign[i] && fatmat[j * (lp->rows + 1) + i] != 0)
printf("% 8.2f ", -fatmat[j * (lp->rows + 1) + i]);
else
printf("% 8.2f ", fatmat[j * (lp->rows + 1) + i]);
if (lp->orig_upbo[i] == lp->infinite)
if (lp->ch_sign[i])
printf(">= ");
else
printf("<= ");
else
printf(" = ");
if (lp->ch_sign[i])
printf("% 8.2f\n", -lp->orig_rh[i]);
else
printf("% 8.2f\n", lp->orig_rh[i]);
}
printf("Type ");
for (i = 1; i <= lp->columns; i++)
if (lp->must_be_int[lp->rows + i] == TRUE)
printf(" Int ");
else
printf(" Real ");
printf("\nupbo ");
for (i = 1; i <= lp->columns; i++)
if (lp->orig_upbo[lp->rows + i] == lp->infinite)
printf(" Inf ");
else
printf("% 8.2f ", lp->orig_upbo[lp->rows + i]);
printf("\nlowbo ");
for (i = 1; i <= lp->columns; i++) printf("% 8.2f ", lp->orig_lowbo[lp->rows + i]);
printf("\n");
for (i = 0; i < lp->nr_lagrange; i++) {
printf("lag[%3d] ", i);
for (j = 1; j <= lp->columns; j++) printf("% 8.2f ", lp->lag_row[i][j]);
if (lp->orig_upbo[i] == lp->infinite)
if (lp->lag_con_type[i] == GE)
printf(">= ");
else if (lp->lag_con_type[i] == LE)
printf("<= ");
else if (lp->lag_con_type[i] == EQ)
printf(" = ");
printf("% 8.2f\n", lp->lag_rhs[i]);
}
free(fatmat);
}
void set_row_name(lprec *lp, int row, nstring new_name) {
int i;
if (!lp->names_used) {
CALLOC(lp->row_name, lp->rows_alloc + 1, nstring);
CALLOC(lp->col_name, lp->columns_alloc + 1, nstring);
lp->names_used = TRUE;
for (i = 0; i <= lp->rows; i++) sprintf(lp->row_name[i], "r_%d", i);
for (i = 1; i <= lp->columns; i++) sprintf(lp->col_name[i], "var_%d", i);
}
strcpy(lp->row_name[row], new_name);
}
void set_col_name(lprec *lp, int column, nstring new_name) {
int i;
if (!lp->names_used) {
CALLOC(lp->row_name, lp->rows_alloc + 1, nstring);
CALLOC(lp->col_name, lp->columns_alloc + 1, nstring);
lp->names_used = TRUE;
for (i = 0; i <= lp->rows; i++) sprintf(lp->row_name[i], "r_%d", i);
for (i = 1; i <= lp->columns; i++) sprintf(lp->col_name[i], "var_%d", i);
}
strcpy(lp->col_name[column], new_name);
}
static REAL minmax_to_scale(REAL min, REAL max) {
REAL scale;
/* should do something sensible when min or max is 0, MB */
if ((min == 0) || (max == 0)) return ((REAL)1);
scale = 1 / pow(10, (log10(min) + log10(max)) / 2);
return (scale);
}
void unscale_columns(lprec *lp) {
int i, j;
/* unscale mat */
for (j = 1; j <= lp->columns; j++)
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++)
lp->mat[i].value /= lp->scale[lp->rows + j];
/* unscale bounds as well */
for (i = lp->rows + 1; i < lp->sum; i++) {
if (lp->orig_lowbo[i] != 0) lp->orig_lowbo[i] *= lp->scale[i];
if (lp->orig_upbo[i] != lp->infinite) lp->orig_upbo[i] *= lp->scale[i];
}
for (i = lp->rows + 1; i <= lp->sum; i++) lp->scale[i] = 1;
lp->columns_scaled = FALSE;
lp->eta_valid = FALSE;
}
void unscale(lprec *lp) {
int i, j;
if (lp->scaling_used) {
/* unscale mat */
for (j = 1; j <= lp->columns; j++)
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++)
lp->mat[i].value /= lp->scale[lp->rows + j];
/* unscale bounds */
for (i = lp->rows + 1; i < lp->sum; i++) {
if (lp->orig_lowbo[i] != 0) lp->orig_lowbo[i] *= lp->scale[i];
if (lp->orig_upbo[i] != lp->infinite) lp->orig_upbo[i] *= lp->scale[i];
}
/* unscale the matrix */
for (j = 1; j <= lp->columns; j++)
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++)
lp->mat[i].value /= lp->scale[lp->mat[i].row_nr];
/* unscale the rhs! */
for (i = 0; i <= lp->rows; i++) lp->orig_rh[i] /= lp->scale[i];
free(lp->scale);
lp->scaling_used = FALSE;
lp->eta_valid = FALSE;
}
}
void auto_scale(lprec *lp) {
int i, j, row_nr, IntUsed;
REAL *row_max, *row_min, *scalechange, absval;
REAL col_max, col_min, col_scale;
if (!lp->scaling_used) {
MALLOC(lp->scale, lp->sum_alloc + 1, REAL);
for (i = 0; i <= lp->sum; i++) lp->scale[i] = 1;
}
MALLOC(row_max, lp->rows + 1, REAL);
MALLOC(row_min, lp->rows + 1, REAL);
MALLOC(scalechange, lp->sum + 1, REAL);
/* initialise min and max values */
for (i = 0; i <= lp->rows; i++) {
row_max[i] = 0;
row_min[i] = lp->infinite;
}
/* calculate min and max absolute values of rows */
for (j = 1; j <= lp->columns; j++) {
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++) {
row_nr = lp->mat[i].row_nr;
absval = my_abs(lp->mat[i].value);
if (absval != 0) {
row_max[row_nr] = my_max(row_max[row_nr], absval);
row_min[row_nr] = my_min(row_min[row_nr], absval);
}
}
}
/* calculate scale factors for rows */
for (i = 0; i <= lp->rows; i++) {
scalechange[i] = minmax_to_scale(row_min[i], row_max[i]);
lp->scale[i] *= scalechange[i];
}
/* now actually scale the matrix */
for (j = 1; j <= lp->columns; j++)
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++)
lp->mat[i].value *= scalechange[lp->mat[i].row_nr];
/* and scale the rhs! */
for (i = 0; i <= lp->rows; i++) lp->orig_rh[i] *= scalechange[i];
free(row_max);
printf("after free r_max. ");
free(row_min);
printf("after free rmin. ");
IntUsed = FALSE;
i = lp->rows + 1;
while (!IntUsed && i <= lp->sum) {
IntUsed = lp->must_be_int[i];
i++;
}
if (!IntUsed) {
REAL col_max, col_min, col_scale;
/* calculate scales */
for (j = 1; j <= lp->columns; j++) {
col_max = 0;
col_min = lp->infinite;
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++) {
if (lp->mat[i].value != 0) {
col_max = my_max(col_max, my_abs(lp->mat[i].value));
col_min = my_min(col_min, my_abs(lp->mat[i].value));
}
}
scalechange[lp->rows + j] = minmax_to_scale(col_min, col_max);
lp->scale[lp->rows + j] *= scalechange[lp->rows + j];
}
/* scale mat */
for (j = 1; j <= lp->columns; j++)
for (i = lp->col_end[j - 1]; i < lp->col_end[j]; i++)
lp->mat[i].value *= scalechange[lp->rows + j];
/* scale bounds as well */
for (i = lp->rows + 1; i < lp->sum; i++) {
if (lp->orig_lowbo[i] != 0) lp->orig_lowbo[i] /= scalechange[i];
if (lp->orig_upbo[i] != lp->infinite) lp->orig_upbo[i] /= scalechange[i];
}
lp->columns_scaled = TRUE;
}
printf("b4 free scch. ");
free(scalechange);
printf("after free scch. ");
lp->scaling_used = TRUE;
lp->eta_valid = FALSE;
}
void reset_basis(lprec *lp) { lp->basis_valid = FALSE; }
void print_solution(lprec *lp) {
int i;
fprintf(stdout, "Value of objective function: %16.10g\n", (double)lp->best_solution[0]);
/* print normal variables */
for (i = 1; i <= lp->columns; i++)
if (lp->names_used)
fprintf(stdout, "%-10s%16.5g\n", lp->col_name[i], (double)lp->best_solution[lp->rows + i]);
else
fprintf(stdout, "Var [%4d] %16.5g\n", i, (double)lp->best_solution[lp->rows + i]);
/* print achieved constraint values */
if (lp->verbose) {
fprintf(stdout, "\nActual values of the constraints:\n");
for (i = 1; i <= lp->rows; i++)
if (lp->names_used)
fprintf(stdout, "%-10s%16.5g\n", lp->row_name[i], (double)lp->best_solution[i]);
else
fprintf(stdout, "%Row [%4d] %16.5g\n", i, (double)lp->best_solution[i]);
}
if ((lp->verbose || lp->print_duals)) {
if (lp->max_level != 1)
fprintf(stdout, "These are the duals from the node that gave the optimal solution.\n");
else
fprintf(stdout, "\nDual values:\n");
for (i = 1; i <= lp->rows; i++)
if (lp->names_used)
fprintf(stdout, "%-10s%16.5g\n", lp->row_name[i], (double)lp->duals[i]);
else
fprintf(stdout, "Row [%4d] %16.5g\n", i, (double)lp->duals[i]);
}
} /* Printsolution */
void write_LP(lprec *lp, FILE *output) {
int i, j;
REAL *row;
nstring var_name;
MALLOC(row, lp->columns + 1, REAL);
if (lp->maximise)
fprintf(output, "max:");
else
fprintf(output, "min:");
get_row(lp, 0, row);
for (i = 1; i <= lp->columns; i++)
if (row[i] != 0) {
if (row[i] == -1)
fprintf(output, " -");
else if (row[i] == 1)
fprintf(output, " +");
else
fprintf(output, " %+g ", row[i]);
if (lp->names_used)
fprintf(output, "%s", lp->col_name[i]);
else
fprintf(output, "x%d", i);
}
fprintf(output, ";\n");
for (j = 1; j <= lp->rows; j++) {
if (lp->names_used) fprintf(output, "%s:", lp->row_name[j]);
get_row(lp, j, row);
for (i = 1; i <= lp->columns; i++)
if (row[i] != 0) {
if (row[i] == -1)
fprintf(output, " -");
else if (row[i] == 1)
fprintf(output, " +");
else
fprintf(output, " %+g ", row[i]);
if (lp->names_used)
fprintf(output, "%s", lp->col_name[i]);
else
fprintf(output, "x%d", i);
}
if (lp->orig_upbo[j] == 0)
fprintf(output, " =");
else if (lp->ch_sign[j])
fprintf(output, " >");
else
fprintf(output, " <");
if (lp->ch_sign[j])
fprintf(output, " %g;\n", -lp->orig_rh[j]);
else
fprintf(output, " %g;\n", lp->orig_rh[j]);
}
for (i = lp->rows + 1; i <= lp->sum; i++) {
if (lp->orig_lowbo[i] != 0) {
if (lp->names_used)
fprintf(output, "%s > %g;\n", lp->col_name[i - lp->rows], lp->orig_lowbo[i]);
else
fprintf(output, "x%d > %g;\n", i - lp->rows, lp->orig_lowbo[i]);
}
if (lp->orig_upbo[i] != lp->infinite) {
if (lp->names_used)
fprintf(output, "%s < %g;\n", lp->col_name[i - lp->rows], lp->orig_upbo[i]);
else
fprintf(output, "x%d < %g;\n", i - lp->rows, lp->orig_upbo[i]);
}
}
i = 1;
while (!lp->must_be_int[lp->rows + i] && i <= lp->columns) i++;
if (i <= lp->columns) {
if (lp->names_used)
fprintf(output, "\nint %s", lp->col_name[i]);
else
fprintf(output, "\nint x%d", i);
i++;
for (; i <= lp->columns; i++)
if (lp->must_be_int[lp->rows + i])
if (lp->names_used)
fprintf(output, ",%s", lp->col_name[i]);
else
fprintf(output, ", x%d", i);
fprintf(output, ";\n");
}
free(row);
}
void write_MPS(lprec *lp, FILE *output) {
int i, j, marker;
REAL *column;
MALLOC(column, lp->rows + 1, REAL);
marker = 0;
fprintf(output, "NAME %s\n", lp->lp_name);
fprintf(output, "ROWS\n");
for (i = 0; i <= lp->rows; i++) {
if (i == 0)
fprintf(output, " N ");
else if (lp->orig_upbo[i] == lp->infinite)
if (lp->ch_sign[i])
fprintf(output, " G ");
else
fprintf(output, " L ");
else
fprintf(output, " E ");
if (lp->names_used)
fprintf(output, "%s\n", lp->row_name[i]);
else
fprintf(output, "r_%d\n", i);
}
fprintf(output, "COLUMNS\n");
j = 0;
for (i = 1; i <= lp->columns; i++) {
if ((lp->must_be_int[i + lp->rows]) && (marker % 2) == 0) {
fprintf(output, " MARK%04d 'MARKER' 'INTORG'\n", marker);
marker++;
}
if ((!lp->must_be_int[i + lp->rows]) && (marker % 2) == 1) {
fprintf(output, " MARK%04d 'MARKER' 'INTEND'\n", marker);
marker++;
}
get_column(lp, i, column);
j = 0;
if (lp->maximise) {
if (column[j] != 0) {
if (lp->names_used)
fprintf(output, " %-8s %-8s %g\n", lp->col_name[i], lp->row_name[j], -column[j]);
else
fprintf(output, " var_%-4d r_%-6d %g\n", i, j, -column[j]);
}
} else {
if (column[j] != 0) {
if (lp->names_used)
fprintf(output, " %-8s %-8s %g\n", lp->col_name[i], lp->row_name[j], column[j]);
else
fprintf(output, " var_%-4d r_%-6d %g\n", i, j, column[j]);
}
}
for (j = 1; j <= lp->rows; j++)
if (column[j] != 0) {
if (lp->names_used)
fprintf(output, " %-8s %-8s %g\n", lp->col_name[i], lp->row_name[j], column[j]);
else
fprintf(output, " var_%-4d r_%-6d %g\n", i, j, column[j]);
}
}
if ((marker % 2) == 1) {
fprintf(output, " MARK%04d 'MARKER' 'INTEND'\n", marker);
marker++;
}
fprintf(output, "RHS\n");
for (i = 1; i <= lp->rows; i++) {
if (lp->ch_sign[i]) {
if (lp->names_used)
fprintf(output, " RHS %-8s %g\n", lp->row_name[i], (double)-lp->orig_rh[i]);
else
fprintf(output, " RHS r_%-6d %g\n", i, (double)-lp->orig_rh[i]);
} else {
if (lp->names_used)
fprintf(output, " RHS %-8s %g\n", lp->row_name[i], (double)lp->orig_rh[i]);
else
fprintf(output, " RHS r_%-6d %g\n", i, (double)lp->orig_rh[i]);
}
}
fprintf(output, "BOUNDS\n");
if (lp->names_used)
for (i = lp->rows + 1; i <= lp->sum; i++) {
if (lp->orig_upbo[i] < lp->infinite)
fprintf(output, " UP BND %-8s %g\n", lp->col_name[i - lp->rows],
(double)lp->orig_upbo[i]);
if (lp->orig_lowbo[i] != 0)
fprintf(output, " LO BND %-8s %g\n", lp->col_name[i - lp->rows],
(double)lp->orig_lowbo[i]);
}
else
for (i = lp->rows + 1; i <= lp->sum; i++) {
if (lp->orig_upbo[i] < lp->infinite)
fprintf(output, " UP BND var_%-4d %g\n", i - lp->rows, (double)lp->orig_upbo[i]);
if (lp->orig_lowbo[i] != 0)
fprintf(output, " LO BND var_%-4d %g\n", i - lp->rows, (double)lp->orig_lowbo[i]);
}
fprintf(output, "ENDATA\n");
free(column);
}
void print_duals(lprec *lp) {
int i;
for (i = 1; i <= lp->rows; i++)
if (lp->names_used)
fprintf(stdout, "%10s [%3d] % 10.4f\n", lp->row_name[i], i, lp->duals[i]);
else
fprintf(stdout, "Dual [%3d] % 10.4f\n", i, lp->duals[i]);
}
void print_scales(lprec *lp) {
int i;
if (lp->scaling_used) {
for (i = 0; i <= lp->rows; i++)
fprintf(stdout, "Row[%3d] scaled at % 10.6f\n", i, lp->scale[i]);
for (i = 1; i <= lp->columns; i++)
fprintf(stdout, "Column[%3d] scaled at % 10.6f\n", i, lp->scale[lp->rows + i]);
}
}
Reference in New Issue View Git Blame Copy Permalink