🚀 go-pugleaf

Submitted-by: Michel Berkelaar <michel@es.ele.tue.nl>
Posting-number: Volume 2, Issue 23
Archive-name: lp_solve/part02

#! /bin/sh
# This is a shell archive.  Remove anything before this line, then unpack
# it by saving it into a file and typing "sh file".  To overwrite existing
# files, type "sh file -c".  You can also feed this as standard input via
# unshar, or by typing "sh <file", e.g..  If this archive is complete, you
# will see the following message at the end:
#		"End of archive 2 (of 2)."
# Contents:  lp_solve/ex5.lp lp_solve/read.c lp_solve/solve.c
# Wrapped by michel@varda.es.ele.tue.nl on Fri Jun 19 09:05:17 1992
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
if test -f lp_solve/ex5.lp -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"lp_solve/ex5.lp\"
else
echo shar: Extracting \"lp_solve/ex5.lp\" \(8873 characters\)
sed "s/^X//" >lp_solve/ex5.lp <<'END_OF_lp_solve/ex5.lp'
X/* correct value of objective function after optimization:  2.3314e+06 */
X
X + 35 COL00001 - 54.9 COL00002 - 54.9 COL00003 + 35 COL00004
X - 23.5 COL00005 - 23.5 COL00006 - 8.72 COL00007 - 9.72 COL00008
X - 6.74 COL00009 - 6.84 COL00010 - 15 COL00011 - 22.5 COL00012 + 500 COL00014
X - 18.7 COL00015 + 258.3 COL00016 + 662 COL00017 - 3 COL00018 - 0.39 COL00019
X - 0.47 COL00020 + 35 COL00021 - 54.9 COL00022 - 54.9 COL00023 + 35 COL00024
X - 23.5 COL00025 - 23.5 COL00026 - 8.72 COL00027 - 9.72 COL00028
X - 6.74 COL00029 - 6.84 COL00030 - 15 COL00031 - 22.5 COL00032 + 500 COL00034
X - 18.7 COL00035 + 258.3 COL00036 + 662 COL00037 - 3 COL00038 - 0.39 COL00039
X - 0.47 COL00040 + 35 COL00041 - 54.9 COL00042 - 54.9 COL00043 + 35 COL00044
X - 23.5 COL00045 - 23.5 COL00046 - 8.72 COL00047 - 9.72 COL00048
X - 6.74 COL00049 - 6.84 COL00050 - 15 COL00051 - 22.5 COL00052 + 500 COL00054
X - 18.7 COL00055 + 258.3 COL00056 + 662 COL00057 - 3 COL00058 - 0.39 COL00059
X - 0.47 COL00060 + 35 COL00061 - 54.9 COL00062 - 54.9 COL00063 + 35 COL00064
X - 23.5 COL00065 - 23.5 COL00066 - 8.72 COL00067 - 9.72 COL00068
X - 6.74 COL00069 - 6.84 COL00070 - 15 COL00071 - 22.5 COL00072 + 500 COL00074
X - 18.7 COL00075 + 258.3 COL00076 + 662 COL00077 - 3 COL00078 - 0.39 COL00079
X - 0.47 COL00080 + 35 COL00081 - 54.9 COL00082 - 54.9 COL00083 + 35 COL00084
X - 23.5 COL00085 - 23.5 COL00086 - 8.72 COL00087 - 9.72 COL00088
X - 6.74 COL00089 - 6.84 COL00090 - 15 COL00091 - 22.5 COL00092 + 500 COL00094
X - 18.7 COL00095 + 258.3 COL00096 + 662 COL00097 - 3 COL00098 - 0.39 COL00099
X - 0.47 COL00100 + 35 COL00101 - 54.9 COL00102 - 54.9 COL00103 + 35 COL00104
X - 23.5 COL00105 - 23.5 COL00106 - 8.72 COL00107 - 9.72 COL00108
X - 6.74 COL00109 - 6.84 COL00110 - 15 COL00111 - 22.5 COL00112 + 500 COL00114
X - 18.7 COL00115 + 258.3 COL00116 + 662 COL00117 - 3 COL00118 - 0.39 COL00119
X - 0.47 COL00120 + 35 COL00121 - 54.9 COL00122 - 54.9 COL00123 + 35 COL00124
X - 23.5 COL00125 - 23.5 COL00126 - 8.72 COL00127 - 9.72 COL00128
X - 6.74 COL00129 - 6.84 COL00130 - 15 COL00131 - 22.5 COL00132 + 500 COL00134
X - 18.7 COL00135 + 258.3 COL00136 + 662 COL00137 - 3 COL00138 - 0.39 COL00139
X - 0.47 COL00140 ;
X COL00001 + COL00002  = 158 ;
X COL00003 + COL00004  = 158 ;
X COL00002 - COL00005  = 0 ;
X COL00003 - COL00006  = 0 ;
X - 0.4 COL00002 - 0.4 COL00003 - 1.7 COL00005 - 1.7 COL00006
X           + 1.5 COL00008 + COL00011  = 3092.96 ;
X COL00008 <= 2566.67 ;
X COL00007 + COL00008 + COL00009 + COL00010 + COL00013  = 6900 ;
X COL00013 <= 1600 ;
X COL00009 + COL00010 >= 800 ;
X COL00005 - COL00014 - COL00015  = 0 ;
X COL00006 - COL00016  = 0 ;
X  - COL00017  = -375.2 ;
X - 0.245 COL00005 - 0.245 COL00006 + 0.2 COL00007 - COL00018
X            = -92.12 ;
X 1.2 COL00009 - COL00019  = -684 ;
X - 0.5 COL00002 - 0.5 COL00003 - 1.4 COL00005 - 1.4 COL00006
X           + COL00010 + COL00012 - COL00020  = -150 ;
X - 0.5 COL00017 + COL00021 + COL00022  = 0 ;
X - 0.5 COL00017 + COL00023 + COL00024  = 0 ;
X COL00022 - COL00025  = 0 ;
X COL00023 - COL00026  = 0 ;
X COL00018 <= 1800 ;
X COL00019 <= 2400 ;
X COL00020 <= 1200 ;
X - 0.48 COL00014 - 0.48 COL00015 - 0.48 COL00016 - 9.32 COL00017
X           - 0.4 COL00022 - 0.4 COL00023 - 1.7 COL00025 - 1.7 COL00026
X           + 1.5 COL00028 + COL00031  = 0 ;
X COL00028 <= 2566.67 ;
X COL00027 + COL00028 + COL00029 + COL00030 + COL00033  = 6900 ;
X COL00033 <= 1600 ;
X COL00029 + COL00030 >= 800 ;
X COL00025 - COL00034 - COL00035  = 0 ;
X COL00026 - COL00036  = 0 ;
X COL00015 + 0.7 COL00017 - COL00037  = 0 ;
X - 0.49 COL00014 - 0.49 COL00015 - 0.49 COL00016 - 0.56 COL00017
X           + COL00018 - 0.245 COL00025 - 0.245 COL00026 + 0.2 COL00027
X           - COL00038  = 0 ;
X - COL00017 + COL00019 + 1.2 COL00029 - COL00039  = 0 ;
X COL00020 - 0.5 COL00022 - 0.5 COL00023 - 1.4 COL00025 - 1.4 COL00026
X           + COL00030 + COL00032 - COL00040  = 0 ;
X COL00017 - COL00037 <= 0 ;
X - 0.5 COL00037 + COL00041 + COL00042  = 0 ;
X - 0.5 COL00037 + COL00043 + COL00044  = 0 ;
X COL00042 - COL00045  = 0 ;
X COL00043 - COL00046  = 0 ;
X COL00038 <= 1800 ;
X COL00039 <= 2400 ;
X COL00040 <= 1200 ;
X - 0.48 COL00034 - 0.48 COL00035 - 0.48 COL00036 - 9.32 COL00037
X           - 0.4 COL00042 - 0.4 COL00043 - 1.7 COL00045 - 1.7 COL00046
X           + 1.5 COL00048 + COL00051  = 0 ;
X COL00048 <= 2566.67 ;
X COL00047 + COL00048 + COL00049 + COL00050 + COL00053  = 6900 ;
X COL00053 <= 1600 ;
X COL00049 + COL00050 >= 800 ;
X COL00045 - COL00054 - COL00055  = 0 ;
X COL00046 - COL00056  = 0 ;
X COL00035 + 0.7 COL00037 - COL00057  = 0 ;
X - 0.49 COL00034 - 0.49 COL00035 - 0.49 COL00036 - 0.56 COL00037
X           + COL00038 - 0.245 COL00045 - 0.245 COL00046 + 0.2 COL00047
X           - COL00058  = 0 ;
X - COL00037 + COL00039 + 1.2 COL00049 - COL00059  = 0 ;
X COL00040 - 0.5 COL00042 - 0.5 COL00043 - 1.4 COL00045 - 1.4 COL00046
X           + COL00050 + COL00052 - COL00060  = 0 ;
X COL00037 - COL00057 <= 0 ;
X - 0.5 COL00057 + COL00061 + COL00062  = 0 ;
X - 0.5 COL00057 + COL00063 + COL00064  = 0 ;
X COL00062 - COL00065  = 0 ;
X COL00063 - COL00066  = 0 ;
X COL00058 <= 1800 ;
X COL00059 <= 2400 ;
X COL00060 <= 1200 ;
X - 0.48 COL00054 - 0.48 COL00055 - 0.48 COL00056 - 9.32 COL00057
X           - 0.4 COL00062 - 0.4 COL00063 - 1.7 COL00065 - 1.7 COL00066
X           + 1.5 COL00068 + COL00071  = 0 ;
X COL00068 <= 2566.67 ;
X COL00067 + COL00068 + COL00069 + COL00070 + COL00073  = 6900 ;
X COL00073 <= 1600 ;
X COL00069 + COL00070 >= 800 ;
X COL00065 - COL00074 - COL00075  = 0 ;
X COL00066 - COL00076  = 0 ;
X COL00055 + 0.7 COL00057 - COL00077  = 0 ;
X - 0.49 COL00054 - 0.49 COL00055 - 0.49 COL00056 - 0.56 COL00057
X           + COL00058 - 0.245 COL00065 - 0.245 COL00066 + 0.2 COL00067
X           - COL00078  = 0 ;
X - COL00057 + COL00059 + 1.2 COL00069 - COL00079  = 0 ;
X COL00060 - 0.5 COL00062 - 0.5 COL00063 - 1.4 COL00065 - 1.4 COL00066
X           + COL00070 + COL00072 - COL00080  = 0 ;
X COL00057 - COL00077 <= 0 ;
X - 0.5 COL00077 + COL00081 + COL00082  = 0 ;
X - 0.5 COL00077 + COL00083 + COL00084  = 0 ;
X COL00082 - COL00085  = 0 ;
X COL00083 - COL00086  = 0 ;
X COL00078 <= 1800 ;
X COL00079 <= 2400 ;
X COL00080 <= 1200 ;
X - 0.48 COL00074 - 0.48 COL00075 - 0.48 COL00076 - 9.32 COL00077
X           - 0.4 COL00082 - 0.4 COL00083 - 1.7 COL00085 - 1.7 COL00086
X           + 1.5 COL00088 + COL00091  = 0 ;
X COL00088 <= 2566.67 ;
X COL00087 + COL00088 + COL00089 + COL00090 + COL00093  = 6900 ;
X COL00093 <= 1600 ;
X COL00089 + COL00090 >= 800 ;
X COL00085 - COL00094 - COL00095  = 0 ;
X COL00086 - COL00096  = 0 ;
X COL00075 + 0.7 COL00077 - COL00097  = 0 ;
X - 0.49 COL00074 - 0.49 COL00075 - 0.49 COL00076 - 0.56 COL00077
X           + COL00078 - 0.245 COL00085 - 0.245 COL00086 + 0.2 COL00087
X           - COL00098  = 0 ;
X - COL00077 + COL00079 + 1.2 COL00089 - COL00099  = 0 ;
X COL00080 - 0.5 COL00082 - 0.5 COL00083 - 1.4 COL00085 - 1.4 COL00086
X           + COL00090 + COL00092 - COL00100  = 0 ;
X COL00077 - COL00097 <= 0 ;
X - 0.5 COL00097 + COL00101 + COL00102  = 0 ;
X - 0.5 COL00097 + COL00103 + COL00104  = 0 ;
X COL00102 - COL00105  = 0 ;
X COL00103 - COL00106  = 0 ;
X COL00098 <= 1800 ;
X COL00099 <= 2400 ;
X COL00100 <= 1200 ;
X - 0.48 COL00094 - 0.48 COL00095 - 0.48 COL00096 - 9.32 COL00097
X           - 0.4 COL00102 - 0.4 COL00103 - 1.7 COL00105 - 1.7 COL00106
X           + 1.5 COL00108 + COL00111  = 0 ;
X COL00108 <= 2566.67 ;
X COL00107 + COL00108 + COL00109 + COL00110 + COL00113  = 6900 ;
X COL00113 <= 1600 ;
X COL00109 + COL00110 >= 800 ;
X COL00105 - COL00114 - COL00115  = 0 ;
X COL00106 - COL00116  = 0 ;
X COL00095 + 0.7 COL00097 - COL00117  = 0 ;
X - 0.49 COL00094 - 0.49 COL00095 - 0.49 COL00096 - 0.56 COL00097
X           + COL00098 - 0.245 COL00105 - 0.245 COL00106 + 0.2 COL00107
X           - COL00118  = 0 ;
X - COL00097 + COL00099 + 1.2 COL00109 - COL00119  = 0 ;
X COL00100 - 0.5 COL00102 - 0.5 COL00103 - 1.4 COL00105 - 1.4 COL00106
X           + COL00110 + COL00112 - COL00120  = 0 ;
X COL00097 - COL00117 <= 0 ;
X - 0.5 COL00117 + COL00121 + COL00122  = 0 ;
X - 0.5 COL00117 + COL00123 + COL00124  = 0 ;
X COL00122 - COL00125  = 0 ;
X COL00123 - COL00126  = 0 ;
X COL00118 <= 1800 ;
X COL00119 <= 2400 ;
X COL00120 <= 1200 ;
X - 0.48 COL00114 - 0.48 COL00115 - 0.48 COL00116 - 9.32 COL00117
X           - 0.4 COL00122 - 0.4 COL00123 - 1.7 COL00125 - 1.7 COL00126
X           + 1.5 COL00128 + COL00131  = 0 ;
X COL00128 <= 2566.67 ;
X COL00127 + COL00128 + COL00129 + COL00130 + COL00133  = 6900 ;
X COL00133 <= 1600 ;
X COL00129 + COL00130 >= 800 ;
X COL00125 - COL00134 - COL00135  = 0 ;
X COL00126 - COL00136  = 0 ;
X COL00115 + 0.7 COL00117 - COL00137  = 0 ;
X - 0.49 COL00114 - 0.49 COL00115 - 0.49 COL00116 - 0.56 COL00117
X           + COL00118 - 0.245 COL00125 - 0.245 COL00126 + 0.2 COL00127
X           - COL00138  = 0 ;
X - COL00117 + COL00119 + 1.2 COL00129 - COL00139  = 0 ;
X COL00120 - 0.5 COL00122 - 0.5 COL00123 - 1.4 COL00125 - 1.4 COL00126
X           + COL00130 + COL00132 - COL00140  = 0 ;
X COL00117 - COL00137 <= 0 ;
X
END_OF_lp_solve/ex5.lp
if test 8873 -ne `wc -c <lp_solve/ex5.lp`; then
    echo shar: \"lp_solve/ex5.lp\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f lp_solve/read.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"lp_solve/read.c\"
else
echo shar: Extracting \"lp_solve/read.c\" \(13928 characters\)
sed "s/^X//" >lp_solve/read.c <<'END_OF_lp_solve/read.c'
X/*
X============================================================================
XNAME    : read.c
XPURPOSE : translation of lp-problem and storage in sparse matrix
XSHORT   : Subroutines for yacc program to store the input in an intermediate
X          data-structure. The yacc and lex programs translate the input.
X          First the problemsize is determined and the date is read into
X          an intermediate structure, then readinput fills the sparse matrix.
XUSAGE   : call yyparse(); to start reading the input.
X          call readinput(); to fill the sparse matrix.
X============================================================================
XRows : contains the amount of rows + 1
X       Rows-1 is the amount of constraints (no bounds)
X       Rows   also contains the rownr 0 which is the objectfunction
XColumns : contains the amount of columns (different variable names
X          found in the constraints)
XNonnuls : contains the amount of nonnuls = sum of different entries
X          of all columns in the constraints and in the objectfunction
XHash_tab : contains all columnnames on the first level of the structure
X           the row information is kept under each column structure
X           in a linked list (also the objext funtion is in this structure)
X           Bound information is also stored under under the column name
XFirst_rside : points to a linked list containing all relational operators
X              and the righthandside values of the constraints
X              the linked list is in reversed order with respect to the
X              rownumbers
X============================================================================
X*/
X#include "defines.h"
X#include "globals.h"
X
X/*---------------------------------------------------------------------------*/
X
X/*
X * errorhandeling routine for yyparse()
X */
Xvoid yyerror(char *string)
X{
X  fprintf(stderr, "PARSING ERROR: %s on line %d, quiting\n", string, yylineno);
X  exit(1);
X}
X
Xvoid check_decl(char *str)
X{
X  if(strcmp(str, "int"))
X    {
X      fprintf(stderr, "Unknown declaration specifier %s on line %d, ignored\n",
X	      str, yylineno);
X      Ignore_decl = TRUE;
X    }
X}
X
Xvoid add_int_var(char *name)
X{
X  hashelem *hp;
X
X  if(Verbose)
X    fprintf(stderr, "int: %s\n", name);
X  if(!(hp = gethash(name)))
X    fprintf(stderr,
X	    "Unknown variable %s declared integer on line %d, ignored\n",
X	    name, yylineno);
X  else
X    hp->must_be_int = 1;
X}
X
X
X/*
X * initialisation of hashtable and globals.
X */
Xvoid init_read(void)
X{
X  int i;
X  Rows = 0;
X  Nonnuls = 0;
X  Columns = 0;
X  for(i = 0; i<HASH_SIZE; Hash_tab[i++] = NULL);
X  CALLOC(First_rside, 1, rside);
X  First_rside->value = 0;
X  First_rside->relat = OF; /* first row (nr 0) is always the objective function */
X} /* init */
X
X/*---------------------------------------------------------------------------*/
X
X/*
X * Returns hashvalue of string
X * hashvalue = sum ( asciivalue of character * indexnumber) mod HASHSIZE
X */
Xint hashvalue(char *string)
X{
X  int i, j;
X  i = j = 0;
X  while((int)string[j]&&j<MAXSTRL)
X    i += (NAMELEN-j)*(int)string[j++];
X  return(i % HASH_SIZE);
X} /* hashvalue */
X
X
X
X/*
X * Returns a pointer to hashelement with colname = variable
X * If this hashelement does not exists, gethash() returns a NULL pointer
X */
Xhashelem *gethash(char *variable) 
X{
X  hashelem *h_tab_p;
X  for(h_tab_p = Hash_tab[hashvalue(variable)];
X      h_tab_p != NULL;
X      h_tab_p = h_tab_p->next)
X    if(strcmp(variable, h_tab_p->colname) == 0)
X      return(h_tab_p);
X  return(h_tab_p);
X} /* gethash */
X
X
X
X/*
X * searchs in column-list (p is pointer to first element of column-list)
X * for column->row = row.
X * getrow() returns a pointer to this column structure.
X * If not found a NULL-pointer is returned
X */
Xcolumn *getrow(column *p,
X	       int row)
X{
X  for(; p != NULL; p = p->next)
X    if(p->row == row)
X      return(p);
X  return(p) ;
X} /* getrow */
X
X
X
X/*
X * Creates a bound record.
X * Set lowbo = 0 and upbo = INFINITE
X *
X */
Xbound *create_bound_rec(void)
X{
X  bound *bp;
X  CALLOC(bp, 1, bound);
X  bp->upbo = INFINITE;
X  bp->lowbo = 0;
X  return(bp);
X} /* create_bound_rec */
X
X
X
X/*
X * clears the tmp_store variable after all information has been copied
X */
Xvoid null_tmp_store(void)
X{
X  tmp_store.value = 0;
X  tmp_store.rhs_value = 0;
X}
X
X/*---------------------------------------------------------------------------*/
X
X/*
X * variable : pointer to text array with name of variable
X * row      : the rownumber of the constraint
X * value    : value of matrixelement
X *            A(row, variable).
X * Sign     : (global)  determines the sign of value.
X * store()  : stores value in matrix
X *	      A(row, variable). If A(row, variable) already contains data,
X *	      value is added to the existing value.
X */
Xvoid store(char *variable,
X	   int row,
X	   double value) 
X{ 
X  hashelem *h_tab_p;
X  column *col_p;
X  int hv;
X
X  if(value == (double)0)
X    return;
X  if((h_tab_p = gethash(variable)) == NULL)
X  {
X    CALLOC(h_tab_p, 1, hashelem);
X    Columns++; /* counter for calloc of final array */
X    hv = hashvalue(variable);
X    h_tab_p->next = Hash_tab[hv];
X    Hash_tab[hv] = h_tab_p;
X    strcpy(h_tab_p->colname, variable);
X    CALLOC(h_tab_p->col, 1, column);
X    Nonnuls++; /* for calloc of final arrays */
X    h_tab_p->col->row = row;
X    h_tab_p->col->value = value;
X  }
X  else
X    if((col_p = getrow(h_tab_p->col, row)) == NULL)
X    {
X      CALLOC(col_p, 1, column);
X      Nonnuls++; /* for calloc of final arrays */
X      col_p->value = value;
X      col_p->row = row;
X      col_p->next = h_tab_p->col;
X      h_tab_p->col = col_p;
X    }
X    else
X      col_p->value += value;
X} /* store */
X
X
X
X/*---------------------------------------------------------------------------*/
X
X/*
X * store relational operator given in yylex[0] in the rightside list.
X * Also checks if it constaint was a bound and if so stores it in the
X * boundslist
X */
Xvoid store_re_op(void)
X{
X  short tmp_relat;
X  switch(yytext[0]) {
X
X  case '=':
X    tmp_relat = EQ;
X    break;
X
X  case '>':
X    tmp_relat=GE;
X    break;
X    
X  case '<':
X    tmp_relat=LE;
X    break;
X    
X  default:
X    break;
X  }
X  
X  if(Lin_term_count > 1) /* it is not a bound */
X    First_rside->relat = tmp_relat;
X  else /* could be a bound */
X    tmp_store.relat = tmp_relat;
X} /* save_re_op */
X
X
X
X/*
X * store RHS value in the rightside structure
X * if type = true then
X */
Xvoid rhs_store(double value)
X{
X  if(Lin_term_count <= 1) /* could be a bound */
X    tmp_store.rhs_value += value;
X  else /* not a bound */
X    First_rside->value += value;
X} /* RHS_store */
X
X
X
X/*
X * store all data in the right place
X * count the amount of lineair terms in a constraint
X * only store in data-structure if the constraint is not a bound
X */
Xvoid var_store(char *var,
X	       int row,
X	       double value)
X{
X  if(strlen(var) > MAXSTRL)
X    {
X      fprintf(stderr, "Variable name too long, at most %d characters allowed",
X	      MAXSTRL);
X      exit(1);
X    }
X  /* also in a bound the same var name can occur more than once. Check for
X     this. Don't increment Lin_term_count */
X
X  if(Lin_term_count != 1 || strcmp(tmp_store.name, var) != 0)
X    Lin_term_count++;
X
X  if(row == 0) /* always store objective function with rownr == 0 */
X    {
X      store(var,  row,  value);
X      return;
X    }
X  
X  if(Lin_term_count == 1) /* don't yet store. could be a bound */
X    {
X      strcpy(tmp_store.name, var);
X      tmp_store.row = row;
X      tmp_store.value += value;
X      return;
X    }
X  
X  if(Lin_term_count == 2) /* now you can also store the first variable */
X    {
X      rside *rp;
X      /* make space for the rhs information */
X      CALLOC(rp, 1, rside);
X      rp->next = First_rside;
X      First_rside = rp;
X      First_rside->value = tmp_store.rhs_value;
X      First_rside->relat = tmp_store.relat;
X      
X      if (tmp_store.value != 0)
X        store(tmp_store.name, tmp_store.row, tmp_store.value);
X
X      null_tmp_store();
X    }
X
X  store(var, row, value);
X} /* var_store */
X
X
X
X/*
X * store the information in tmp_store because it is a bound
X */
Xstore_bounds(void)
X{
X  if (tmp_store.value != 0)
X    {
X      hashelem *h_tab_p;
X      int hv;
X
X      if((h_tab_p = gethash(tmp_store.name)) == NULL)
X        {
X          /* a new columnname is found, create an entry in the hashlist */
X          CALLOC(h_tab_p, 1, hashelem);
X          Columns++; /* counter for calloc of final array */
X	  hv            = hashvalue(tmp_store.name);
X	  h_tab_p->next = Hash_tab[hv];
X          Hash_tab[hv]  = h_tab_p;
X          strcpy(h_tab_p->colname, tmp_store.name);
X          /* create a place to store bounds information */
X          h_tab_p->bnd = create_bound_rec();
X        }
X      else
X        if(h_tab_p->bnd == NULL)
X          /* create a place to store bounds information */
X          h_tab_p->bnd = create_bound_rec();
X      /* else bound_rec already exists */
X   
X      if(tmp_store.value < 0) /* divide by negative number, */
X  	    	              /* relational operator may change */
X        {
X          if (tmp_store.relat == GE)
X            tmp_store.relat = LE;
X          else if (tmp_store.relat == LE)
X            tmp_store.relat = GE;
X        }
X 
X      if ((tmp_store.relat == GE) || (tmp_store.relat == EQ))
X        h_tab_p->bnd->lowbo = tmp_store.rhs_value / tmp_store.value;
X      if ((tmp_store.relat == LE) || (tmp_store.relat == EQ))
X        h_tab_p->bnd->upbo = tmp_store.rhs_value / tmp_store.value;
X    }
X  
X  null_tmp_store();
X} /* store_bounds */
X
X/* ========================================================================= */
X
X
X/*
X * reallocates eta
X */
Xvoid resize_eta(void)
X{
X  Cur_eta_size *= 2;
X  if (Verbose)
X    printf("Resizing Eta_value and Eta_rownr, new size is %d\n", Cur_eta_size);
X  if(!(Eta_value = realloc(Eta_value, Cur_eta_size * sizeof(double))))
X    {
X      fprintf(stderr, "Error: cannot realloc Eta_value to size %d (entries)\n",
X	      Cur_eta_size);
X      exit(1);
X    }
X  if(!(Eta_rownr = realloc(Eta_rownr, Cur_eta_size * sizeof(int))))
X    {
X      fprintf(stderr, "Error: cannot realloc Eta_rownr to size %d (entries)\n",
X	      Cur_eta_size);
X      exit(1);
X    }
X} /* resize_eta */
X
X
X
X
X/*
X * transport the data from the intermediate structure to the sparse matrix
X * and free the intermediate structure
X */
Xvoid readinput(int *cend,
X	       double *rh,
X	       short *relat,
X	       double *lowbo,
X	       double *upbo,
X	       matrec *mat,
X	       nstring *names)
X{
X  int    i, j, k, index, nn_ind;
X  column *cp,*tcp;    /* tcp (temporary cp) points to memory-space to free */
X  hashelem *hp,*thp;
X  bound *bp;
X  int   x;
X  rside *rp;
X  intrec *irp;
X  
X  /* initialize lower and upper bound arrays */
X  for (i = 0; i <= Sum; i++)
X    {
X      lowbo[i] = 0;
X      upbo[i] = INFINITE;
X    }
X  
X  /* fill names with the rownames */
X  for (i = 0; i <= Rows; i++)
X    {
X      /* the first row is row zero (the objective function)             */
X      sprintf(names[i], "%s%d", STD_ROW_NAME_PREFIX, i);
X    }
X  
X  for (i = Rows;i >= 0;i--)
X    {
X      rp = First_rside;
X      relat[i] = rp->relat;
X      rh[i] = rp->value;
X      First_rside = rp->next;
X      free(rp);  /* free memory when data has been read */
X    }
X  
X  /* change upperbound to zero if the relational operator is the equal sign */
X  for (i = 1; i <= Rows; i++)
X    if (relat[i] == EQ)
X      upbo[i] = 0;
X  
X  /* start reading the Hash_list structure */
X  index = 0;
X  nn_ind = 0;
X  
X  for (i = 0;i < HASH_SIZE;i++)
X    {
X      hp = Hash_tab[i];
X      while (hp != NULL)
X        {
X	  /* put an index in the cend array when a new name is found */
X	  cend[index++] = nn_ind;
X
X	  /* check if it must be an integer variable */
X	  if(hp->must_be_int)
X	    {
X	      CALLOC(irp, 1, intrec);
X	      irp->varnr = Rows + index;
X	      irp->next  = First_int;
X	      First_int  = irp;
X	    }
X	  /* check for bound */
X	  if (hp->bnd != NULL)
X	    {
X	      bp = hp->bnd;
X              lowbo[Rows+index] = bp->lowbo;
X              upbo[Rows+index] = bp->upbo;
X              free(bp); /* free memory when data has been read*/
X	    }
X	  
X	  /* copy name of column variable */
X	  strcpy(names[Rows+index], hp->colname);
X
X	  /* put matrix values in sparse matrix */
X          cp = hp->col;
X          while (cp!=NULL)
X            {
X              mat[nn_ind].rownr = cp->row;
X	      mat[nn_ind].value = cp->value;
X	      nn_ind++;
X	      tcp = cp;
X	      cp = cp->next;
X	      free(tcp);  /* free memory when data has been read */
X	    }
X	  thp = hp;
X	  hp = hp->next;
X          free(thp);    /* free memory when data has been read */
X        }
X    }
X  cend[index] = nn_ind; 
X  
X  if (Verbose)
X    {
X      printf("\n");
X      printf("**********Data read**********\n");
X      printf("Rows    : %d\n", Rows);
X      printf("Columns : %d\n", Columns);
X      printf("Nonnuls : %d\n", Nonnuls);
X      printf("\nSparse Matrix\nRow right hand side\n");
X      printf("%4s  %8s %3s %s\n\n", "nr", "row", "rel", "Value");
X      for (i = 0; i <= Rows; i++)
X	{
X	  printf("%4d  %8s ", i, names[i]);
X	  if (relat[i] == LE) printf(" < ");
X	  if (relat[i] == EQ) printf(" = ");
X	  if (relat[i] == GE) printf(" > ");
X	  if (relat[i] == OF) printf("ObF");
X	  printf(" %f\n", rh[i]);
X	}
X      
X      printf("\nMatrix contents\n%8s  %8s  %s\n\n","colname","row","value");
X      j = 0;
X      for (i = 0; i < Nonnuls; i++)
X	{
X	  if (i == cend[j])
X	    printf("%8s  %8s  %f\n", names[Rows+ ++j], names[mat[i].rownr],
X		   mat[i].value);
X	  else
X	    printf("          %8s  %e\n", names[mat[i].rownr], mat[i].value);
X	}
X      printf("\nBounds\n%8s %3s %s\n\n", "name", "rel", "value");
X      for (i = 0; i <= Sum; i++)
X	{
X	  if (upbo[i] < (INFINITE/10000))  /* double is to small to contain
X					      INFINITE thus is rounded */
X	    printf("%8s  <  %f\n", names[i], upbo[i]);
X	  if (lowbo[i] > 0)
X	    printf("%8s  >  %f\n", names[i], lowbo[i]);
X	}
X      
X      printf("\n**********End data**********\n\n");
X    }
X} /* readinput */
X
X
X/* ===================== END OF read.c ===================================== */
X
END_OF_lp_solve/read.c
if test 13928 -ne `wc -c <lp_solve/read.c`; then
    echo shar: \"lp_solve/read.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
if test -f lp_solve/solve.c -a "${1}" != "-c" ; then 
  echo shar: Will not over-write existing file \"lp_solve/solve.c\"
else
echo shar: Extracting \"lp_solve/solve.c\" \(18876 characters\)
sed "s/^X//" >lp_solve/solve.c <<'END_OF_lp_solve/solve.c'
X#include "defines.h"
X#include "globals.h"
X
X
Xvoid  condensecol(int rownr,
X		  int numc,
X		  double *pcol)
X{
X  int  i, elnr;
X  
X  if (Verbose)
X    printf("condensecol\n");
X  elnr = Endetacol[numc];
X
X  if (elnr + Rows + 2 > Cur_eta_size) /* maximum local growth of Eta */
X    resize_eta();
X
X  for (i = 0; i <= Rows; i++)
X    if (i != rownr && pcol[i] != 0)
X      {
X	Eta_rownr[elnr] = i;
X	Eta_value[elnr] = pcol[i];
X	elnr++;
X      }
X  Eta_rownr[elnr] = rownr;
X  Eta_value[elnr] = pcol[rownr];
X  elnr++;
X  Endetacol[numc + 1] = elnr;
X} /* condensecol */
X
Xvoid  addetacol(int *numc)
X{
X  int    i, j, k;
X  double theta;
X  
X  if (Verbose)
X    printf("addetacol\n");
X  j = Endetacol[(*numc)];
X  (*numc)++;
X  k = Endetacol[(*numc)];
X  theta = 1 / (double) Eta_value[k - 1];
X  Eta_value[k - 1] = theta;
X  for (i = j; i < k - 1; i++)
X    Eta_value[i] *= -theta;
X} /* addetacol */
X
X
Xvoid  presolve(void)
X{
X  int    i, j, rownr;
X  double theta;
X  int    *num, *rownum;
X
X  CALLOC(num, Rows + 2, int);
X  CALLOC(rownum, Rows + 2, int);
X  
X  if (Verbose)
X    printf("presolve\n");
X  Rh[0] = 0;
X  for (i = 0; i <= Sum; i++)
X    Lower[i] = TRUE;
X  for (i = 1; i <= Rows; i++)
X    Basis[i] = TRUE;
X  for (i = 1; i <= Columns; i++)
X    Basis[Rows + i] = FALSE;
X  for (i = 0; i <= Rows; i++)
X    Bas[i] = i;
X
X  /* we want to maximise the objective function, so: */
X  Chsign[0] = TRUE;
X
X  for (i = 1; i <= Rows; i++)
X    if ((Relat[i] == GE) && (Upbo[i] == INFINITE))
X      Chsign[i] = TRUE;
X    else
X      Chsign[i] = FALSE;
X
X  /* invert values in rows with Chsign[i] == TRUE */
X  for (i = 1; i <= Rows; i++)
X    if (Chsign[i])
X      Rh[i] = -Rh[i];
X  for (i = 0; i < Nonnuls; i++)
X    if (Chsign[Mat[i].rownr])
X      Mat[i].value = -Mat[i].value;
X
X  /* handle lower bounds by transforming to 0 */
X  for (i = 1; i <= Columns; i++)
X    if (Lowbo[Rows + i] > 0)
X      {
X	theta = Lowbo[Rows + i];
X	Upbo[Rows + i] -= theta;
X	for (j = Cend[i - 1]; j < Cend[i]; j++)
X	  Rh[Mat[j].rownr] -= theta * Mat[j].value;
X      }
X
X  for (i = 1; i <= Rows; i++)
X    {
X      num[i] = 0;
X      rownum[i] = 0;
X    }
X  for (i = 0; i < Nonnuls; i++)
X    rownum[Mat[i].rownr]++;
X  Rend[0] = 0;
X  for (i = 1; i <= Rows; i++)
X    Rend[i] = Rend[i - 1] + rownum[i];
X  for (i = 1; i <= Columns; i++)
X    for (j = Cend[i - 1]; j < Cend[i]; j++)
X      {
X	rownr = Mat[j].rownr;
X	if (rownr != 0)
X	  {
X	    num[rownr]++;
X	    Colno[Rend[rownr - 1] + num[rownr]] = i;
X	  }
X      }
X  free((char *)num);
X  free((char *)rownum);
X} /* presolve */
X
Xvoid  setpivcol(short *lower, 
X		int varin,
X		int numeta,
X		double *pcol)
X{
X  int    i, colnr;
X  double f;
X  
X  if (Verbose)
X    printf("setpivcol\n");
X  if ((*lower))
X    f = 1;
X  else
X    f = -1;
X  for (i = 0; i <= Rows; i++)
X    pcol[i] = 0;
X  if (varin > Rows)
X    {
X      colnr = varin - Rows;
X      for (i = Cend[colnr - 1]; i < Cend[colnr]; i++)
X	pcol[Mat[i].rownr] = Mat[i].value * f;
X      pcol[0] -= Extrad * f;
X    }
X  else
X    if ((*lower))
X      pcol[varin] = 1;
X    else
X      pcol[varin] = -1;
X  ftran(1, numeta, pcol);
X} /* setpivcol */
X
X
Xshort  colprim(int *colnr,
X	       int *numeta,
X	       short *minit,
X	       double *drow)
X{
X  int    varnr, i, j;
X  double f, dpiv;
X  
X  if (Verbose)
X    printf("colprim\n");
X  dpiv = -EPSD;
X  (*colnr) = 0;
X  if (!(*minit))
X    {
X      for (i = 1; i <= Rows; i++)
X	drow[i] = 0;
X      drow[0] = 1;
X      btran((*numeta), drow);
X      for (i = 1; i <= Columns; i++)
X	{
X	  varnr = Rows + i;
X	  if (!Basis[varnr])
X	    if (Upbo[varnr] > 0)
X	      {
X		f = 0;
X		for (j = Cend[i - 1]; j < Cend[i]; j++)
X		  f += drow[Mat[j].rownr] * Mat[j].value;
X		drow[varnr] = f;
X	      }
X	}
X      for (i = 1; i <= Sum; i++)
X	if (abs(drow[i]) < EPSD)
X	  drow[i] = 0;
X    }
X  for (i = 1; i <= Sum; i++)
X    if (!Basis[i])
X      if (Upbo[i] > 0)
X	{
X	  if (Lower[i])
X	    f = drow[i];
X	  else
X	    f = -drow[i];
X	  if (f < dpiv)
X	    {
X	      dpiv = f;
X	      (*colnr) = i;
X	    }
X	}
X  return ((*colnr) > 0);
X} /* colprim */
X
X
Xvoid  minoriteration(int colnr,
X		     int rownr,
X		     int *numeta)
X{
X  int    i, j, k, wk, varin, varout, elnr;
X  double pivot, theta;
X  
X  if (Verbose)
X    printf("minoriteration\n");
X  varin = colnr + Rows;
X  elnr = Endetacol[(*numeta)];
X  wk = elnr;
X  (*numeta)++;
X  if (Extrad != 0)
X    {
X      Eta_rownr[elnr] = 0;
X      Eta_value[elnr] = -Extrad;
X      elnr++;
X    }
X  for (j = Cend[colnr - 1] ; j < Cend[colnr]; j++)
X    {
X      k = Mat[j].rownr;
X      if (k == 0 && Extrad != 0)
X	Eta_value[Endetacol[(*numeta) - 1]] += Mat[j].value;
X      else
X	if (k != rownr)
X	  {
X	    Eta_rownr[elnr] = k;
X	    Eta_value[elnr] = Mat[j].value;
X	    elnr++;
X	  }
X	else
X	  pivot = Mat[j].value;
X    }
X  Eta_rownr[elnr] = rownr;
X  Eta_value[elnr] = 1 / (double) pivot;
X  elnr++;
X  theta = Rhs[rownr] / (double) pivot;
X  Rhs[rownr] = theta;
X  for (i = wk; i < elnr - 1; i++)
X    Rhs[Eta_rownr[i]] -= theta * Eta_value[i];
X  varout = Bas[rownr];
X  Bas[rownr] = varin;
X  Basis[varout] = FALSE;
X  Basis[varin] = TRUE;
X  for (i = wk; i < elnr - 1; i++)
X    Eta_value[i] /= - (double) pivot;
X  /* printf("minoriteration: new etaindex: %d\n", elnr); */
X  Endetacol[(*numeta)] = elnr;
X} /* minoriteration */
X
X
Xvoid  rhsmincol(double *theta,
X		int rownr,
X		int varin,
X		int numeta)
X{
X  int    i, j, k, varout;
X  double f;
X  
X  if (Verbose)
X    printf("rhsmincol\n");
X  if (rownr > Rows + 1)
X    {
X      fprintf(stderr, "Error: rhsmincol called with rownr: %d, Rows: %d\n",
X	      rownr, Rows);
X    }
X  j = Endetacol[numeta];
X  k = Endetacol[numeta + 1];
X  for (i = j; i < k; i++)
X    {
X      f = Rhs[Eta_rownr[i]] - (*theta) * Eta_value[i];
X      if (abs(f) < EPSB)
X	Rhs[Eta_rownr[i]] = 0;
X      else
X	Rhs[Eta_rownr[i]] = f;
X    }
X  Rhs[rownr] = (*theta);
X  varout = Bas[rownr];
X  Bas[rownr] = varin;
X  Basis[varout] = FALSE;
X  Basis[varin] = TRUE;
X} /* rhsmincol */
X
Xvoid  invert(int n,
X	     int *numeta,
X	     int *numinv)
X{
X  int    i, j, v, wk, numit, varnr, rownr, colnr, varin;
X  double f, theta;
X  double *pcol;
X  short  *frow;
X  short  *fcol;
X  int    *rownum, *col, *row;
X  int    *colnum;
X  short  *lbas;
X  double *rhs1;
X 
X  CALLOC(rownum, Rows + 1, int);
X  CALLOC(col, Rows + 1, int);
X  CALLOC(row, Rows + 1, int);
X  CALLOC(pcol, Rows + 2, double);
X  CALLOC(frow, Rows + 2, short);
X  CALLOC(fcol, Columns + 1, short);
X  CALLOC(colnum, Columns + 2, int);
X  CALLOC(lbas, Sum +2, short);
X  CALLOC(rhs1, Rows + 2, double);
X 
X  if (Verbose)
X    printf("a: invert:%d %d  %12f\n", (*numeta), Endetacol[(*numeta)], Rhs[0]);
X  
X  for(i = 0; i <= Sum + 1; i++)
X    lbas[i] = Basis[i];
X  for(i = 0; i <= Rows + 1; i++)
X    rhs1[i] = Rhs[i];
X  for (i = 0; i <= Rows; i++)
X    frow[i] = TRUE;
X  for (i = 0; i < n; i++)
X    fcol[i] = FALSE;
X  for (i = 0; i < Rows; i++)
X    rownum[i] = 0;
X  for (i = 1; i <= n; i++)
X    colnum[i] = 0;
X  for (i = 0; i <= Rows; i++)
X    if (Bas[i] > Rows)
X      fcol[Bas[i] - Rows - 1] = TRUE;
X    else
X      frow[Bas[i]] = FALSE;
X  for (i = 1; i <= Rows; i++)
X    if (frow[i])
X      for (j = Rend[i - 1] + 1; j <= Rend[i]; j++)
X	{
X	  wk = Colno[j];
X	  if (fcol[wk - 1])
X	    {
X	      colnum[wk]++;
X	      rownum[i - 1]++;
X	    }
X	}
X  for (i = 1; i <= Rows; i++)
X    Bas[i] = i;
X  for (i = 1; i <= Rows; i++)
X    Basis[i] = TRUE;
X  for (i = 1; i <= n; i++)
X    Basis[i + Rows] = FALSE;
X  for (i = 0; i <= Rows; i++)
X    Rhs[i] = Rh[i];
X  for (i = 1; i <= n; i++)
X    {
X      varnr = Rows + i;
X      if (!Lower[varnr])
X	{
X	  theta = Upbo[varnr];
X	  for (j = Cend[i - 1]; j < Cend[i]; j++)
X	    Rhs[Mat[j].rownr] -= theta * Mat[j].value;
X	}
X    }
X  for (i = 1; i <= Rows; i++)
X    if (!Lower[i])
X      Rhs[i] -= Upbo[i];
X  (*numeta) = 0;
X  v = 0;
X  rownr = 0;
X  (*numinv) = 0;
X  numit = 0;
X  while (v < Rows)
X    {
X      rownr++;
X      if (rownr > Rows)
X	rownr = 1;
X      v++;
X      if (rownum[rownr - 1] == 1)
X	if (frow[rownr])
X	  {
X	    v = 0;
X	    j = Rend[rownr - 1] + 1;
X	    while (!(fcol[Colno[j] - 1]))
X	      j++;
X	    colnr = Colno[j];
X	    fcol[colnr - 1] = FALSE;
X	    colnum[colnr] = 0;
X	    for (j = Cend[colnr - 1]; j < Cend[colnr]; j++)
X	      if (frow[Mat[j].rownr])
X		rownum[Mat[j].rownr - 1]--;
X	    frow[rownr] = FALSE;
X	    minoriteration(colnr, rownr, numeta);
X	  }
X    }
X  v = 0;
X  colnr = 0;
X  while (v < n)
X    {
X      colnr++;
X      if (colnr > n)
X	colnr = 1;
X      v++;
X      if (colnum[colnr] == 1)
X	if (fcol[colnr - 1])
X	  {
X	    v = 0;
X	    j = Cend[colnr - 1] + 1;
X	    while (!(frow[Mat[j - 1].rownr]))
X	      j++;
X	    rownr = Mat[j - 1].rownr;
X	    frow[rownr] = FALSE;
X	    rownum[rownr - 1] = 0;
X	    for (j = Rend[rownr - 1] + 1; j <= Rend[rownr]; j++)
X	      if (fcol[Colno[j] - 1])
X		colnum[Colno[j]]--;
X	    fcol[colnr - 1] = FALSE;
X	    numit++;
X	    col[numit - 1] = colnr;
X	    row[numit - 1] = rownr;
X	  }
X    }
X  for (j = 1; j <= n; j++)
X    if (fcol[j - 1])
X      {
X	fcol[j - 1] = FALSE;
X	setpivcol(&Lower[Rows + j], j + Rows, (*numeta), pcol);
X	rownr = 1;
X	while (!(frow[rownr] && pcol[rownr]))
X	  rownr++; /* this sometimes sets rownr to Rows + 2 and makes */
X		   /* rhsmincol crash. MB */
X	frow[rownr] = FALSE;
X	condensecol(rownr, (*numeta), pcol);
X	theta = Rhs[rownr] / (double) pcol[rownr];
X	rhsmincol(&theta, rownr, Rows + j, (*numeta));
X	addetacol(numeta);
X      }
X  for (i = numit - 1; i >= 0; i--)
X    {
X      colnr = col[i];
X      rownr = row[i];
X      varin = colnr + Rows;
X      for (j = 0; j <= Rows; j++)
X	pcol[j] = 0;
X      for (j = Cend[colnr - 1]; j < Cend[colnr]; j++)
X	pcol[Mat[j].rownr] = Mat[j].value;
X      pcol[0] -= Extrad;
X      condensecol(rownr, (*numeta), pcol);
X      theta = Rhs[rownr] / (double) pcol[rownr];
X      rhsmincol(&theta, rownr, varin, (*numeta));
X      addetacol(numeta);
X    }
X  for (i = 1; i <= Rows; i++)
X    if (abs(Rhs[i]) < EPSB)
X      Rhs[i] = 0;
X  f = 0;
X  for (i = 1; i <= Rows; i++)
X    if (Rhs[i] < 0)
X      f += Rhs[i];
X    else
X      if (Rhs[i] > Upbo[Bas[i]])
X        f = f + Upbo[Bas[i]] - Rhs[i];
X  if (Verbose)
X    printf("b: invert:%d %d  %12f  %12f\n", (*numeta), Endetacol[(*numeta)],
X	   Rhs[0], f);
X
X  free((char *)rownum);
X  free((char *)col);
X  free((char *)row);
X  free((char *)pcol);
X  free((char *)frow);
X  free((char *)fcol);
X  free((char *)colnum);
X  free((char *)lbas);
X  free((char *)rhs1);
X} /* invert */
X
X
Xshort  rowprim(int *rownr,
X	       double *theta,
X	       double *pcol)
X{
X  int    i;
X  double f, quot;
X  
X  if (Verbose)
X    printf("rowprim\n");
X  (*rownr) = 0;
X  (*theta) = INFINITE;
X  for (i = 1; i <= Rows; i++)
X    {
X      f = pcol[i];
X      if (f != 0)
X	{
X	  quot = (*theta) * 10;
X	  if (f > 0)
X	    quot = Rhs[i] / (double) f;
X	  else
X	    if (Upbo[Bas[i]] < INFINITE)
X	      quot = (Rhs[i] - Upbo[Bas[i]]) / (double) f;
X	  if (quot < (*theta))
X	    {
X	      (*theta) = quot;
X	      (*rownr) = i;
X	    }
X	}
X    }
X  return ((*rownr) > 0);
X} /* rowprim */
X
X
X
Xvoid  iteration(int *numc,
X		int *rownr,
X		int *varin,
X		double *theta,
X		double *up,
X		short *minit,
X		short *low,
X		short *primair,
X		short *smotes,
X		int *iter,
X		int *numinv)
X{
X  int    i, k, varout;
X  double f, pivot;
X  
X  if (Verbose)
X    printf("iteration\n");
X  (*iter)++;
X  (*minit) = (*theta) > (*up);
X  if ((*minit))
X    {
X      (*theta) = (*up);
X      (*low) = !(*low);
X    }
X  k = Endetacol[(*numc) + 1];
X  pivot = Eta_value[k - 1];
X  for (i = Endetacol[(*numc)]; i < k; i++)
X    {
X      f = Rhs[Eta_rownr[i]] - (*theta) * Eta_value[i];
X      if (abs(f) < EPSB)
X	Rhs[Eta_rownr[i]] = 0;
X      else
X	Rhs[Eta_rownr[i]] = f;
X    }
X  if (!(*minit))
X    {
X      Rhs[(*rownr)] = (*theta);
X      varout = Bas[(*rownr)];
X      Bas[(*rownr)] = (*varin);
X      Basis[varout] = FALSE;
X      Basis[(*varin)] = TRUE;
X      if ((*primair) && pivot < 0)
X	Lower[varout] = FALSE;
X      if (!(*low) && (*up) != INFINITE)
X	{
X	  (*low) = TRUE;
X	  Rhs[(*rownr)] = (*up) - Rhs[(*rownr)];
X	  for (i = Endetacol[(*numc)]; i < k; i++)
X	    Eta_value[i] = -Eta_value[i];
X	}
X      addetacol(numc);
X      (*numinv)++;
X    }
X  if ((*smotes) && Verbose)
X    {
X      printf("Iteration %d: ", (*iter));
X      if ((*minit))
X	{
X	  printf("%4d", (*varin));
X	  if (Lower[(*varin)])
X	    printf("  ");
X	  else
X	    printf(" u");
X	  printf("%22c", ' ');
X	}
X      else
X	{
X	  printf("%4d", (*varin));
X	  if (Lower[(*varin)])
X	    printf("  ");
X	  else
X	    printf(" u");
X	  printf("%4d", varout);
X	  if (Lower[varout])
X	    printf("  ");
X	  else
X	    printf(" u");
X	  printf("%12.5f", pivot);
X	}
X      printf("    %12f", Pcol[0]);
X      if (!(*primair))
X	{
X	  f = 0;
X	  for (i = 1; i <= Rows; i++)
X	    if (Rhs[i] < 0)
X	      f += Rhs[i];
X	    else
X	      if (Rhs[i] > Upbo[Bas[i]])
X		f = f + Upbo[Bas[i]] - Rhs[i];
X	  printf("  %12f", f);
X	}
X      else
X	printf("    %12f", Rhs[0]);
X      printf("\n");
X    }
X} /* iteration */
X
X
Xint  solvelp(void)
X{
X  int    i, j, iter, varnr;
X  double f, theta;
X  short  primair, doiter, doinvert, smotes;
X  double *drow, *prow;
X  int    numinv, invnum;
X  int    numeta;
X  short  artif, minit;
X  int    colnr, rownr;
X 
X  CALLOC(drow, Sum + 2, double);
X  CALLOC(prow, Sum + 2, double);
X 
X  numinv = 0;
X  invnum = 50; /* number of iterations between inversions */
X  numeta = 0;
X  iter = 0;
X  minit = FALSE;
X  primair = TRUE;
X  i = 0;
X  smotes = TRUE;
X
X  while (i != Rows && primair)
X    {
X      i++;
X      primair = Rhs[i] >= 0 && Rhs[i] <= Upbo[Bas[i]];
X    }
X  
X  if (!primair)
X    {
X      drow[0] = 1;
X      for (i = 1; i <= Rows; i++)
X	drow[i] = 0;
X      Extrad = 0;
X      for (i = 1; i <= Columns; i++)
X	{
X	  varnr = Rows + i;
X	  drow[varnr] = 0;
X	  for (j = Cend[i - 1]; j < Cend[i]; j++)
X	    if (drow[Mat[j].rownr] != 0)
X	      drow[varnr] += drow[Mat[j].rownr] * Mat[j].value;
X	  if (drow[varnr] < Extrad)
X	    Extrad = drow[varnr];
X	}
X      artif = Extrad < -EPSD;
X    }
X  else
X    {
X      artif = FALSE;
X      Extrad = 0;
X    }
X  if (Verbose)
X    printf("artificial:%12.5f\n", Extrad);
X  minit = FALSE;
X  do {
X    doiter = FALSE;
X    doinvert = FALSE;
X    if (primair)
X      {
X	rownr = 0;
X	if (colprim(&colnr, &numeta, &minit, drow))
X	  {
X	    setpivcol(&Lower[colnr], colnr, numeta, Pcol);
X	    if (rowprim(&rownr, &theta, Pcol))
X	      {
X		doiter = TRUE;
X		condensecol(rownr, numeta, Pcol);
X	      }
X	  }
X      }
X    else
X      {
X	if (!minit)
X	  rowdual(&rownr);
X	if (!doinvert)
X	  if (rownr > 0)
X	    {
X	      if (coldual(&numeta, &rownr, &colnr, &minit, prow, drow))
X		{
X		  doiter = TRUE;
X		  setpivcol(&Lower[colnr], colnr, numeta, Pcol);
X		  condensecol(rownr, numeta, Pcol);
X		  f = Rhs[rownr] - Upbo[Bas[rownr]];
X		  if (f > 0)
X		    {
X		      theta = f / (double) Pcol[rownr];
X		      if (theta <= Upbo[colnr])
X			Lower[Bas[rownr]] = !Lower[Bas[rownr]];
X		    }
X		  else /* getting div by zero here ... MB */
X		    theta = Rhs[rownr] / (double) Pcol[rownr];
X		}
X	    }
X	  else
X	    if (artif)
X	      {
X		primair = TRUE;
X		artif = FALSE;
X		printf("primal\n");
X		doinvert = TRUE;
X		Extrad = 0;
X	      }
X	    else
X	      colnr = 0;
X      }
X    if (doiter)
X      iteration(&numeta, &rownr, &colnr, &theta, &Upbo[colnr],
X		&minit, &Lower[colnr], &primair, &smotes, &iter, &numinv);
X    if (numinv >= invnum)
X      doinvert = TRUE;
X    if (doinvert)
X      invert(Columns, &numeta, &numinv);
X  } while (rownr && colnr || doinvert);
X
X  free((char *)drow);
X  free((char *)prow);
X
X  return(rownr || colnr); /* return nonzero if anything was wrong */
X} /* solvelp */
X
X
Xvoid construct_solution(double *sol)
X{
X  int    i, j, basi;
X  double f;
X
X  /* zero all results of rows */
X  bzero((char *)sol, (Rows + 1) * sizeof(double));
X
X  for (i = Rows + 1; i <= Sum; i++)
X    sol[i] = Lowbo[i];
X  for (i = 1; i <= Rows; i++)
X    {
X      basi = Bas[i];
X      if (basi > Rows)
X	sol[basi] += Rhs[i];
X    }
X  for (i = Rows + 1; i <= Sum; i++)
X    if (!Basis[i] && !Lower[i])
X      sol[i] += Upbo[i];
X  for (j = 1; j <= Columns; j++)
X    {
X      f = sol[Rows + j];
X      if (f != 0)
X	for (i = Cend[j - 1]; i < Cend[j]; i++)
X	  sol[Mat[i].rownr] += f * Mat[i].value;
X    }
X  
X  for (i = 0; i <= Rows; i++)
X    {
X      if (abs(sol[i]) < EPSB)
X	sol[i] = 0;
X      else
X	if (Chsign[i])
X	  sol[i] = -sol[i];
X    }
X} /* construct_solution */
X
X
Xint solve(double *upbo,
X	  double *lowbo)
X{
X  int    i, failure, notint;
X  intrec *ptr;
X
X  Level++;
X  debug_print("starting solve");
X
X  /* make fresh copies of Mat, Rh, Upbo, Lowbo, as solving changes them */
X  bcopy((char*)Orig_mat,  (char *)Mat,   (Nonnuls + 1) * sizeof(matrec));
X  bcopy((char*)Orig_rh,   (char *)Rh,    (Rows + 2)    * sizeof(double));
X  bcopy((char*)upbo,      (char *)Upbo,  (Sum + 2)     * sizeof(double));
X  bcopy((char*)lowbo,     (char *)Lowbo, (Sum + 2)     * sizeof(double));
X
X  Rend[0] = 0;
X  presolve();
X
X  for (i = 0; i <= Rows; i++)
X    Rhs[i] = Rh[i];
X
X  Endetacol[0] = 0;
X  failure = solvelp();
X
X  if(failure)
X    debug_print("this problem has no solution");
X
X  if(!failure) /* there is a solution */
X    {
X      construct_solution(Solution);
X
X      debug_print("a solution was found");
X      debug_print_solution(Solution);
X
X      /* if this solution is worse than the best sofar, this branch must die */
X      if(Solution[0] <= Best_solution[0])
X	{
X	  debug_print("but it was worse than the best sofar, discarded");
X	  Level--;
X	  return(TRUE);
X	}
X
X      /* check if solution contains enough ints */
X      notint = 0;
X      for(ptr = First_int; !notint && ptr; ptr = ptr->next)
X	if(!is_int(Solution[ptr->varnr]))
X	  notint = ptr->varnr;
X
X      if(notint) /* there is at least one value not yet int */
X	{
X	  /* set up two new problems */
X	  double *new_upbo, *new_lowbo;
X	  double new_bound;
X	  int    res1, res2;
X
X	  /* allocate room for them */
X	  MALLOC(new_upbo,  Sum + 2, double);
X	  MALLOC(new_lowbo, Sum + 2, double);
X
X	  bcopy((char *)upbo,  (char *)new_upbo,  (Sum + 2) * sizeof(double));
X	  bcopy((char *)lowbo, (char *)new_lowbo, (Sum + 2) * sizeof(double));
X
X	  debug_print("not enough ints. Selecting var %s, val: %10.3g",
X		      Names[notint], Solution[notint]);
X
X	  new_bound = floor(Solution[notint]);
X	  new_upbo[notint] = new_bound;
X
X	  debug_print("starting first subproblem with bounds:");
X	  debug_print_bounds(new_upbo, lowbo);
X
X	  res1 = solve(new_upbo, lowbo);
X
X	  new_bound += 1;
X	  new_lowbo[notint] = new_bound;
X
X	  debug_print("starting second subproblem with bounds:");
X	  debug_print_bounds(upbo, new_lowbo);
X
X	  res2 = solve(upbo, new_lowbo);
X
X	  if(res1 && res2) /* both failed */
X	    failure = TRUE;
X	  else
X	    failure = FALSE;
X
X	  FREE(new_upbo);
X	  FREE(new_lowbo);
X	}
X      else /* all required values are int */
X	{
X	  debug_print("--> Valid solution found");
X	  if(Solution[0] > Best_solution[0]) /* Current solution better */
X	    {
X	      debug_print("the best sofar. Prev: %10.3g, New: %10.3g",
X			  Best_solution[0], Solution[0]);
X	      bcopy((char *)Solution, (char *)Best_solution,
X		    (Sum + 1) * sizeof(double));
X	      if(Show_results)
X		{
X		  fprintf(stderr, "Intermediate solution:\n");
X		  print_solution(stderr, Best_solution);
X		}
X	    }
X	}
X    }
X  Level--;
X  return(failure);
X} /* solve */
END_OF_lp_solve/solve.c
if test 18876 -ne `wc -c <lp_solve/solve.c`; then
    echo shar: \"lp_solve/solve.c\" unpacked with wrong size!
fi
# end of overwriting check
fi
echo shar: End of archive 2 \(of 2\).
cp /dev/null ark2isdone
MISSING=""
for I in 1 2 ; do
    if test ! -f ark${I}isdone ; then
	MISSING="${MISSING} ${I}"
    fi
done
if test "${MISSING}" = "" ; then
    echo You have unpacked both archives.
    rm -f ark[1-9]isdone
else
    echo You still need to unpack the following archives:
    echo "        " ${MISSING}
fi
##  End of shell archive.
exit 0

exit 0 # Just in case...
-- 
        Andrew Patrick acting as Comp.Sources.Reviewed Moderator
              Department of Communications, Ottawa, CANADA
                         csr@calvin.dgbt.doc.CA