LibRCG  3.1.1
rlp.c
Go to the documentation of this file.
1 
9 #include <stdio.h>
10 #include <float.h>
11 #include <math.h>
12 #include "rlp.h"
13 
18 #define POS(R,C,NC) ((R)*(NC)+(C))
19 
31 static void fmprint(double* matrix,int nrows,int ncols,FILE* file)
32 {
33  int i,j;
34  fprintf(file,"\n +");
35  for(i=0;i<ncols;i++)
36  {
37  fprintf(file,"--------------+");
38  }
39  fputs("\n",file);
40  for(i=0;i<nrows;i++)
41  {
42  fprintf(file," |");
43  for(j=0;j<ncols;j++)
44  {
45  fprintf(file," %12.4f |",matrix[POS(i,j,ncols)]);
46  }
47  fputs("\n",file);
48  }
49  fprintf(file," +");
50  for(i=0;i<ncols;i++)
51  {
52  fprintf(file,"--------------+");
53  }
54  fputs("\n\n",file);
55 }
56 
57 //==============================================================================
58 
72 static double minimumc(double* matrix,int nrows,int ncols,int col,int* row)
73 {
74  int i;
75  double res=DBL_MAX;
76  for(i=0;i<nrows;i++)
77  {
78  if(matrix[POS(i,col,ncols)]<res)
79  {
80  res=matrix[POS(i,col,ncols)];
81  *row=i;
82  }
83  }
84  return res;
85 }
86 
87 //==============================================================================
88 
101 static double minimumr(double* matrix,int ncols,int row,int* col)
102 {
103  int i;
104  double res=DBL_MAX;
105  for(i=0;i<ncols;i++)
106  {
107  if(matrix[POS(row,i,ncols)]<res)
108  {
109  res=matrix[POS(row,i,ncols)];
110  *col=i;
111  }
112  }
113  return res;
114 }
115 
116 //==============================================================================
117 
118 int simplex(double* a,int n,int m,FILE* file)
119 {
120  double aux;
121  int pos=0,nm2,err=0;
122  nm2=n+m+2;
123  aux=minimumc(&a[nm2],m,nm2,n+m,&pos);
124  if(file) fmprint(a,m+1,nm2,file);
125  if(aux<0) err=simplexd(a,n,m,pos+1,file);
126  if(!err)
127  {
128  aux=minimumr(a,nm2,0,&pos);
129  if(aux<0) err=simplexp(a,n,m,pos,file);
130  }
131  return err;
132 }
133 
134 //==============================================================================
135 
136 int simplexp(double* a,int n,int m,int pos,FILE* file)
137 {
138  int nm,pivotc,pivotr,i,j,k,result=0;
139  double pivot,coef,aux,r;
140  nm=n+m;
141  k=nm+2;
142  aux=-1;
143  pivotc=pos;
144  while(aux<0)
145  {
146  pivotr=-1;
147  aux=DBL_MAX;
148  for(i=1;i<m+1;i++)
149  {
150  if(a[POS(i,pivotc,k)]>0)
151  {
152  r=a[POS(i,nm,k)]/a[POS(i,pivotc,k)];
153  if(r<aux)
154  {
155  aux=r;
156  pivotr=i;
157  }
158  }
159  }
160  if(pivotr==-1) result=1;
161  else
162  {
163  pivot=a[POS(pivotr,pivotc,k)];
164  for(i=0;i<nm+1;i++)
165  {
166  a[POS(pivotr,i,k)]/=pivot;
167  }
168  for(i=1;i<m+1;i++)
169  {
170  if(i!=pivotr)
171  {
172  coef=-a[POS(i,pivotc,k)];
173  for(j=0;j<nm+1;j++)
174  a[POS(i,j,k)]+=coef*a[POS(pivotr,j,k)];
175  }
176  }
177  aux=DBL_MAX;
178  coef=-a[POS(0,pivotc,k)];
179  a[POS(pivotr,nm+1,k)]=pivotc;
180  for(i=0;i<nm+1;i++)
181  {
182  a[POS(0,i,k)]+=coef*a[POS(pivotr,i,k)];
183  if(a[POS(0,i,k)]<aux)
184  {
185  aux=a[POS(0,i,k)];
186  pivotc=i;
187  }
188  }
189  if(file) fmprint(a,m+1,k,file);
190  }
191  }
192  return result;
193 }
194 
195 //==============================================================================
196 
197 int simplexd(double* a,int n,int m,int pos,FILE* file)
198 {
199  int pivotc,pivotr,nm,i,j,k,result=0;
200  double aux,pivot,r,coef;
201  nm=n+m;
202  k=nm+2;
203  aux=-1;
204  pivotr=pos;
205  while(aux<0)
206  {
207  pivotc=-1;
208  aux=DBL_MAX;
209 
210  for(i=0;i<nm;i++)
211  {
212  if(a[POS(pivotr,i,k)]<0)
213  {
214  r=a[POS(0,i,k)]/a[POS(pivotr,i,k)];
215  if(fabs(r)<aux)
216  {
217  aux=r;
218  pivotc=i;
219  }
220  }
221  }
222  if(pivotc==-1) result=1;
223  else
224  {
225  pivot=a[POS(pivotr,pivotc,k)];
226  for(i=0;i<nm+1;i++)
227  {
228  a[POS(pivotr,i,k)]/=pivot;
229  }
230  aux=a[POS(pivotr,nm,k)];
231  for(i=0;i<m+1;i++)
232  {
233  if(i!=pivotr)
234  {
235  coef=a[POS(i,pivotc,k)]>0?a[POS(i,pivotc,k)]:-a[POS(i,pivotc,k)];
236  for(j=0;j<nm+1;j++)
237  {
238  a[POS(i,j,k)]+=coef*a[POS(pivotr,j,k)];
239  }
240  if(a[POS(i,nm,k)]<aux)
241  {
242  aux=a[POS(i,nm,k)];
243  pivotr=i;
244  }
245  }
246  }
247  a[POS(pivotr,nm+1,k)]=pivotc;
248  if(file) fmprint(a,m+1,k,file);
249  }
250  }
251  return result;
252 }
Implementation of linear programming functions.
int simplexp(double *a, int n, int m, int pos, FILE *file)
Applies the Simplex algorithm to a primal optimization problem.
Definition: rlp.c:136
#define POS(R, C, NC)
Given a row (R), a column (C), and the number of columns (NC) of a matrix, computes an equivalent 1D ...
Definition: rlp.c:18
static double minimumc(double *matrix, int nrows, int ncols, int col, int *row)
Finds the minimum value of a matrix's column.
Definition: rlp.c:72
static void fmprint(double *matrix, int nrows, int ncols, FILE *file)
Prints a matrix associated to an optimization problem.
Definition: rlp.c:31
static double minimumr(double *matrix, int ncols, int row, int *col)
Finds the minimum value of a matrix's row.
Definition: rlp.c:101
int simplex(double *a, int n, int m, FILE *file)
Applies the Simplex Algorithm to an optimization problem.
Definition: rlp.c:118
int simplexd(double *a, int n, int m, int pos, FILE *file)
Applies the Simplex algorithm to a dual optimization problem.
Definition: rlp.c:197

LibRCG © 2004-2015   Rui Carlos Gonçalves