1
2 /*
3 Mixed integer programming optimization drivers for lp_solve v5.0+
4 ----------------------------------------------------------------------------------
5 Author: Michel Berkelaar (to lp_solve v3.2),
6 Kjell Eikland
7 Contact:
8 License terms: LGPL.
9
10 Requires: stdarg.h, lp_lib.h
11
12 Release notes:
13 v5.0.0 3 1 January 2004 New unit isolating reporting routines.
14 v5.2.0.0 1 December 2005 Addition of Matrix Market writing function.
15
16 ----------------------------------------------------------------------------------
17 */
18
19 #include <stdio.h>
20 #include <stdarg.h>
21 #include <string.h>
22
23 #include "lp_lib.h"
24 #include "lp_scale.h"
25 #include "commonlib.h"
26 #include "lp_report.h"
27
28 #include <R.h>
29 #include <R_ext/Print.h>
30
31
32 #include "mmio.h"
33
34 #ifdef FORTIFY
35 # include "lp_fortify.h"
36 #endif
37
38 /* Define buffer-size controled function mapping */
39 # if defined _MSC_VER
40 # define vsnprintf _vsnprintf
41 # endif
42
43 /* Various reporting functions for lp_solve */
44 /* ------------------------------------------------------------------------- */
45
46 /* First define general utilties for reporting and output */
explain(lprec * lp,char * format,...)47 char * __VACALL explain(lprec *lp, char *format, ...)
48 {
49 char buff[DEF_STRBUFSIZE+1];
50 va_list ap;
51
52 va_start(ap, format);
53 vsnprintf(buff, DEF_STRBUFSIZE, format, ap);
54 allocCHAR(lp, &(lp->ex_status), (int) strlen(buff), AUTOMATIC);
55 strcpy(lp->ex_status, buff);
56 va_end(ap);
57 return( lp->ex_status );
58 }
report(lprec * lp,int level,char * format,...)59 void __VACALL report(lprec *lp, int level, char *format, ...)
60 {
61 static char buff[DEF_STRBUFSIZE+1];
62 static va_list ap;
63
64 if(lp == NULL) {
65 va_start(ap, format);
66 REvprintf( format, ap);
67 va_end(ap);
68 }
69 else if(level <= lp->verbose) {
70 va_start(ap, format);
71 if(lp->writelog != NULL) {
72 vsnprintf(buff, DEF_STRBUFSIZE, format, ap);
73 lp->writelog(lp, lp->loghandle, buff);
74 }
75 if(lp->outstream != NULL) {
76 vfprintf(lp->outstream, format, ap);
77 /* if(lp->outstream != stdout) */
78 fflush(lp->outstream);
79 }
80 va_end(ap);
81 }
82 #ifdef xParanoia
83 if(level == CRITICAL)
84 raise(SIGSEGV);
85 #endif
86 }
87
print_indent(lprec * lp)88 STATIC void print_indent(lprec *lp)
89 {
90 int i;
91
92 report(lp, NEUTRAL, "%2d", lp->bb_level);
93 if(lp->bb_level < 50) /* useless otherwise */
94 for(i = lp->bb_level; i > 0; i--)
95 report(lp, NEUTRAL, "--");
96 else
97 report(lp, NEUTRAL, " *** too deep ***");
98 report(lp, NEUTRAL, "> ");
99 } /* print_indent */
100
debug_print(lprec * lp,char * format,...)101 STATIC void debug_print(lprec *lp, char *format, ...)
102 {
103 va_list ap;
104
105 if(lp->bb_trace) {
106 print_indent(lp);
107 va_start(ap, format);
108 if (lp == NULL)
109 {
110 REvprintf( format, ap);
111 /* fputc('\n', stderr); */
112 }
113 else if(lp->debuginfo != NULL)
114 {
115 char buff[DEF_STRBUFSIZE+1];
116 vsnprintf(buff, DEF_STRBUFSIZE, format, ap);
117 lp->debuginfo(lp, lp->loghandle, buff);
118 }
119 va_end(ap);
120 }
121 } /* debug_print */
122
debug_print_solution(lprec * lp)123 STATIC void debug_print_solution(lprec *lp)
124 {
125 int i;
126
127 if(lp->bb_trace)
128 for (i = lp->rows + 1; i <= lp->sum; i++) {
129 print_indent(lp);
130 report(lp, NEUTRAL, "%s " RESULTVALUEMASK "\n",
131 get_col_name(lp, i - lp->rows),
132 (double)lp->solution[i]);
133 }
134 } /* debug_print_solution */
135
debug_print_bounds(lprec * lp,REAL * upbo,REAL * lowbo)136 STATIC void debug_print_bounds(lprec *lp, REAL *upbo, REAL *lowbo)
137 {
138 int i;
139
140 if(lp->bb_trace)
141 for(i = lp->rows + 1; i <= lp->sum; i++) {
142 if(lowbo[i] == upbo[i]) {
143 print_indent(lp);
144 report(lp, NEUTRAL, "%s = " RESULTVALUEMASK "\n", get_col_name(lp, i - lp->rows),
145 (double)lowbo[i]);
146 }
147 else {
148 if(lowbo[i] != 0) {
149 print_indent(lp);
150 report(lp, NEUTRAL, "%s > " RESULTVALUEMASK "\n", get_col_name(lp, i - lp->rows),
151 (double)lowbo[i]);
152 }
153 if(upbo[i] != lp->infinite) {
154 print_indent(lp);
155 report(lp, NEUTRAL, "%s < " RESULTVALUEMASK "\n", get_col_name(lp, i - lp->rows),
156 (double)upbo[i]);
157 }
158 }
159 }
160 } /* debug_print_bounds */
161
162 /* List a vector of LREAL values for the given index range */
blockWriteLREAL(FILE * output,char * label,LREAL * vector,int first,int last)163 void blockWriteLREAL(FILE *output, char *label, LREAL *vector, int first, int last)
164 {
165 int i, k = 0;
166
167 fprintf(output, "%s", label);
168 fprintf(output, "\n");
169 for(i = first; i <= last; i++) {
170 fprintf(output, " %18g", vector[i]);
171 k++;
172 if(my_mod(k, 4) == 0) {
173 fprintf(output, "\n");
174 k = 0;
175 }
176 }
177 if(my_mod(k, 4) != 0)
178 fprintf(output, "\n");
179 }
180
181 /* List the current user data matrix columns over the selected row range */
blockWriteAMAT(FILE * output,const char * label,lprec * lp,int first,int last)182 void blockWriteAMAT(FILE *output, const char *label, lprec* lp, int first, int last)
183 {
184 int i, j, k = 0;
185 int nzb, nze, jb;
186 double hold;
187 MATrec *mat = lp->matA;
188
189 if(!mat_validate(mat))
190 return;
191 if(first < 0)
192 first = 0;
193 if(last < 0)
194 last = lp->rows;
195
196 fprintf(output, "%s", label);
197 fprintf(output, "\n");
198
199 if(first == 0) {
200 for(j = 1; j <= lp->columns; j++) {
201 hold = get_mat(lp, 0, j);
202 fprintf(output, " %18g", hold);
203 k++;
204 if(my_mod(k, 4) == 0) {
205 fprintf(output, "\n");
206 k = 0;
207 }
208 }
209 if(my_mod(k, 4) != 0) {
210 fprintf(output, "\n");
211 k = 0;
212 }
213 first++;
214 }
215 nze = mat->row_end[first-1];
216 for(i = first; i <= last; i++) {
217 nzb = nze;
218 nze = mat->row_end[i];
219 if(nzb >= nze)
220 jb = lp->columns+1;
221 else
222 jb = ROW_MAT_COLNR(nzb);
223 for(j = 1; j <= lp->columns; j++) {
224 if(j < jb)
225 hold = 0;
226 else {
227 hold = get_mat(lp, i, j);
228 nzb++;
229 if(nzb < nze)
230 jb = ROW_MAT_COLNR(nzb);
231 else
232 jb = lp->columns+1;
233 }
234 fprintf(output, " %18g", hold);
235 k++;
236 if(my_mod(k, 4) == 0) {
237 fprintf(output, "\n");
238 k = 0;
239 }
240 }
241 if(my_mod(k, 4) != 0) {
242 fprintf(output, "\n");
243 k = 0;
244 }
245 }
246 if(my_mod(k, 4) != 0)
247 fprintf(output, "\n");
248 }
249
250 /* List the current basis matrix columns over the selected row range */
blockWriteBMAT(FILE * output,const char * label,lprec * lp,int first,int last)251 void blockWriteBMAT(FILE *output, const char *label, lprec* lp, int first, int last)
252 {
253 int i, j, jb, k = 0;
254 double hold;
255
256 if(first < 0)
257 first = 0;
258 if(last < 0)
259 last = lp->rows;
260
261 fprintf(output, "%s", label);
262 fprintf(output, "\n");
263
264 for(i = first; i <= last; i++) {
265 for(j = 1; j <= lp->rows; j++) {
266 jb = lp->var_basic[j];
267 if(jb <= lp->rows) {
268 if(jb == i)
269 hold = 1;
270 else
271 hold = 0;
272 }
273 else
274 hold = get_mat(lp, i, j);
275 if(i == 0)
276 modifyOF1(lp, jb, &hold, 1);
277 hold = unscaled_mat(lp, hold, i, jb);
278 fprintf(output, " %18g", hold);
279 k++;
280 if(my_mod(k, 4) == 0) {
281 fprintf(output, "\n");
282 k = 0;
283 }
284 }
285 if(my_mod(k, 4) != 0) {
286 fprintf(output, "\n");
287 k = 0;
288 }
289 }
290 if(my_mod(k, 4) != 0)
291 fprintf(output, "\n");
292 }
293
294 /* Do a generic readable data dump of key lp_solve model variables;
295 principally for run difference and debugging purposes */
REPORT_debugdump(lprec * lp,char * filename,MYBOOL livedata)296 MYBOOL REPORT_debugdump(lprec *lp, char *filename, MYBOOL livedata)
297 {
298 /* FILE *output = stdout; */
299 FILE *output;
300 MYBOOL ok;
301
302 ok = (MYBOOL) ((filename == NULL) || ((output = fopen(filename,"w")) != NULL));
303 if(!ok)
304 return(ok);
305 if((filename == NULL) && (lp->outstream != NULL))
306 output = lp->outstream;
307
308 fprintf(output, "\nGENERAL INFORMATION\n-------------------\n\n");
309 fprintf(output, "Model size: %d rows (%d equalities, %d Lagrangean), %d columns (%d integers, %d SC, %d SOS, %d GUB)\n",
310 lp->rows, lp->equalities, get_Lrows(lp), lp->columns,
311 lp->int_vars, lp->sc_vars, SOS_count(lp), GUB_count(lp));
312 fprintf(output, "Data size: %d model non-zeros, %d invB non-zeros (engine is %s)\n",
313 get_nonzeros(lp), my_if(lp->invB == NULL, 0, lp->bfp_nonzeros(lp, FALSE)), lp->bfp_name());
314 fprintf(output, "Internal sizes: %d rows allocated, %d columns allocated, %d columns used, %d eta length\n",
315 lp->rows_alloc, lp->columns_alloc, lp->columns, my_if(lp->invB == NULL, 0, lp->bfp_colcount(lp)));
316 fprintf(output, "Memory use: %d sparse matrix, %d eta\n",
317 lp->matA->mat_alloc, my_if(lp->invB == NULL, 0, lp->bfp_memallocated(lp)));
318 fprintf(output, "Parameters: Maximize=%d, Names used=%d, Scalingmode=%d, Presolve=%d, SimplexPivot=%d\n",
319 is_maxim(lp), lp->names_used, lp->scalemode, lp->do_presolve, lp->piv_strategy);
320 fprintf(output, "Precision: EpsValue=%g, EpsPrimal=%g, EpsDual=%g, EpsPivot=%g, EpsPerturb=%g\n",
321 lp->epsvalue, lp->epsprimal, lp->epsdual, lp->epspivot, lp->epsperturb);
322 fprintf(output, "Stability: AntiDegen=%d, Improvement=%d, Split variables at=%g\n",
323 lp->improve, lp->anti_degen, lp->negrange);
324 fprintf(output, "B&B settings: BB pivot rule=%d, BB branching=%s, BB strategy=%d, Integer precision=%g, MIP gaps=%g,%g\n",
325 lp->bb_rule, my_boolstr(lp->bb_varbranch), lp->bb_floorfirst, lp->epsint, lp->mip_absgap, lp->mip_relgap);
326
327 fprintf(output, "\nCORE DATA\n---------\n\n");
328 blockWriteINT(output, "Column starts", lp->matA->col_end, 0, lp->columns);
329 blockWriteINT(output, "row_type", lp->row_type, 0, lp->rows);
330 blockWriteREAL(output, "orig_rhs", lp->orig_rhs, 0, lp->rows);
331 blockWriteREAL(output, "orig_lowbo", lp->orig_lowbo, 0, lp->sum);
332 blockWriteREAL(output, "orig_upbo", lp->orig_upbo, 0, lp->sum);
333 blockWriteINT(output, "row_type", lp->row_type, 0, lp->rows);
334 blockWriteBOOL(output, "var_type", lp->var_type, 0, lp->columns, TRUE);
335 blockWriteAMAT(output, "A", lp, 0, lp->rows);
336
337 if(livedata) {
338 fprintf(output, "\nPROCESS DATA\n------------\n\n");
339 blockWriteREAL(output, "Active rhs", lp->rhs, 0, lp->rows);
340 blockWriteINT(output, "Basic variables", lp->var_basic, 0, lp->rows);
341 blockWriteBOOL(output, "is_basic", lp->is_basic, 0, lp->sum, TRUE);
342 blockWriteREAL(output, "lowbo", lp->lowbo, 0, lp->sum);
343 blockWriteREAL(output, "upbo", lp->upbo, 0, lp->sum);
344 if(lp->scalars != NULL)
345 blockWriteREAL(output, "scalars", lp->scalars, 0, lp->sum);
346 }
347
348 if(filename != NULL)
349 fclose(output);
350 return(ok);
351 }
352
353
354 /* High level reports for model results */
355
REPORT_objective(lprec * lp)356 void REPORT_objective(lprec *lp)
357 {
358 if(lp->outstream == NULL)
359 return;
360 fprintf(lp->outstream, "\nValue of objective function: %g\n",
361 (double)lp->best_solution[0]);
362 fflush(lp->outstream);
363 }
364
REPORT_solution(lprec * lp,int columns)365 void REPORT_solution(lprec *lp, int columns)
366 {
367 int i, j, n;
368 REAL value;
369 presolveundorec *psundo = lp->presolve_undo;
370 MYBOOL NZonly = (MYBOOL) ((lp->print_sol & AUTOMATIC) > 0);
371
372 if(lp->outstream == NULL)
373 return;
374
375 fprintf(lp->outstream, "\nActual values of the variables:\n");
376 if(columns <= 0)
377 columns = 2;
378 n = 0;
379 for(i = 1; i <= psundo->orig_columns; i++) {
380 j = psundo->orig_rows + i;
381 value = get_var_primalresult(lp, j);
382 if(NZonly && (fabs(value) < lp->epsprimal))
383 continue;
384 n = (n+1) % columns;
385 fprintf(lp->outstream, "%-20s %12g", get_origcol_name(lp, i), (double) value);
386 if(n == 0)
387 fprintf(lp->outstream, "\n");
388 else
389 fprintf(lp->outstream, " ");
390 }
391
392 fflush(lp->outstream);
393 } /* REPORT_solution */
394
REPORT_constraints(lprec * lp,int columns)395 void REPORT_constraints(lprec *lp, int columns)
396 {
397 int i, n;
398 REAL value;
399 MYBOOL NZonly = (MYBOOL) ((lp->print_sol & AUTOMATIC) > 0);
400
401 if(lp->outstream == NULL)
402 return;
403
404 if(columns <= 0)
405 columns = 2;
406
407 fprintf(lp->outstream, "\nActual values of the constraints:\n");
408 n = 0;
409 for(i = 1; i <= lp->rows; i++) {
410 value = (double)lp->best_solution[i];
411 if(NZonly && (fabs(value) < lp->epsprimal))
412 continue;
413 n = (n+1) % columns;
414 fprintf(lp->outstream, "%-20s %12g", get_row_name(lp, i), value);
415 if(n == 0)
416 fprintf(lp->outstream, "\n");
417 else
418 fprintf(lp->outstream, " ");
419 }
420
421 fflush(lp->outstream);
422 }
423
REPORT_duals(lprec * lp)424 void REPORT_duals(lprec *lp)
425 {
426 int i;
427 REAL *duals, *dualsfrom, *dualstill, *objfrom, *objtill, *objfromvalue;
428 MYBOOL ret;
429
430 if(lp->outstream == NULL)
431 return;
432
433 ret = get_ptr_sensitivity_objex(lp, &objfrom, &objtill, &objfromvalue, NULL);
434 if(ret) {
435 fprintf(lp->outstream, "\nObjective function limits:\n");
436 fprintf(lp->outstream, " From Till FromValue\n");
437 for(i = 1; i <= lp->columns; i++)
438 if(!is_splitvar(lp, i))
439 fprintf(lp->outstream, "%-20s %15.7g %15.7g %15.7g\n", get_col_name(lp, i),
440 (double)objfrom[i - 1], (double)objtill[i - 1], (double)objfromvalue[i - 1]);
441 }
442
443 ret = get_ptr_sensitivity_rhs(lp, &duals, &dualsfrom, &dualstill);
444 if(ret) {
445 fprintf(lp->outstream, "\nDual values with from - till limits:\n");
446 fprintf(lp->outstream, " Dual value From Till\n");
447 for(i = 1; i <= lp->sum; i++)
448 fprintf(lp->outstream, "%-20s %15.7g %15.7g %15.7g\n",
449 (i <= lp->rows) ? get_row_name(lp, i) : get_col_name(lp, i - lp->rows),
450 (double)duals[i - 1], (double)dualsfrom[i - 1], (double)dualstill[i - 1]);
451 fflush(lp->outstream);
452 }
453 }
454
455 /* Printing of sensitivity analysis reports */
REPORT_extended(lprec * lp)456 void REPORT_extended(lprec *lp)
457 {
458 int i, j;
459 REAL hold;
460 REAL *duals, *dualsfrom, *dualstill, *objfrom, *objtill;
461 MYBOOL ret;
462
463 ret = get_ptr_sensitivity_obj(lp, &objfrom, &objtill);
464 report(lp, NORMAL, " \n");
465 report(lp, NORMAL, "Primal objective:\n");
466 report(lp, NORMAL, " \n");
467 report(lp, NORMAL, " Column name Value Objective Min Max\n");
468 report(lp, NORMAL, " --------------------------------------------------------------------------\n");
469 for(j = 1; j <= lp->columns; j++) {
470 hold = get_mat(lp,0,j);
471 report(lp, NORMAL, " %-25s " MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK "\n",
472 get_col_name(lp,j),
473 my_precision(hold,lp->epsprimal),
474 my_precision(hold*lp->best_solution[lp->rows+j],lp->epsprimal),
475 my_precision((ret) ? objfrom[j - 1] : 0.0,lp->epsprimal),
476 my_precision((ret) ? objtill[j - 1] : 0.0,lp->epsprimal));
477 }
478 report(lp, NORMAL, " \n");
479
480 ret = get_ptr_sensitivity_rhs(lp, &duals, &dualsfrom, &dualstill);
481 report(lp, NORMAL, "Primal variables:\n");
482 report(lp, NORMAL, " \n");
483 report(lp, NORMAL, " Column name Value Slack Min Max\n");
484 report(lp, NORMAL, " --------------------------------------------------------------------------\n");
485 for(j = 1; j <= lp->columns; j++)
486 report(lp, NORMAL, " %-25s " MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK "\n",
487 get_col_name(lp,j),
488 my_precision(lp->best_solution[lp->rows+j],lp->epsprimal),
489 my_precision(my_inflimit(lp, (ret) ? duals[lp->rows+j-1] : 0.0),lp->epsprimal),
490 my_precision((ret) ? dualsfrom[lp->rows+j-1] : 0.0,lp->epsprimal),
491 my_precision((ret) ? dualstill[lp->rows+j-1] : 0.0,lp->epsprimal));
492
493 report(lp, NORMAL, " \n");
494 report(lp, NORMAL, "Dual variables:\n");
495 report(lp, NORMAL, " \n");
496 report(lp, NORMAL, " Row name Value Slack Min Max\n");
497 report(lp, NORMAL, " --------------------------------------------------------------------------\n");
498 for(i = 1; i <= lp->rows; i++)
499 report(lp, NORMAL, " %-25s " MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK MPSVALUEMASK "\n",
500 get_row_name(lp,i),
501 my_precision((ret) ? duals[i - 1] : 0.0, lp->epsprimal),
502 my_precision(lp->best_solution[i], lp->epsprimal),
503 my_precision((ret) ? dualsfrom[i - 1] : 0.0,lp->epsprimal),
504 my_precision((ret) ? dualstill[i - 1] : 0.0,lp->epsprimal));
505
506 report(lp, NORMAL, " \n");
507 }
508
509 /* A more readable lp-format report of the model; antiquated and not updated */
REPORT_lp(lprec * lp)510 void REPORT_lp(lprec *lp)
511 {
512 int i, j;
513
514 if(lp->outstream == NULL)
515 return;
516
517 if(lp->matA->is_roworder) {
518 report(lp, IMPORTANT, "REPORT_lp: Cannot print lp while in row entry mode.\n");
519 return;
520 }
521
522 fprintf(lp->outstream, "Model name: %s\n", get_lp_name(lp));
523 fprintf(lp->outstream, " ");
524
525 for(j = 1; j <= lp->columns; j++)
526 fprintf(lp->outstream, "%8s ", get_col_name(lp,j));
527
528 fprintf(lp->outstream, "\n%simize ", (is_maxim(lp) ? "Max" : "Min"));
529 for(j = 1; j <= lp->columns; j++)
530 fprintf(lp->outstream, "%8g ", get_mat(lp, 0, j));
531 fprintf(lp->outstream, "\n");
532
533 for(i = 1; i <= lp->rows; i++) {
534 fprintf(lp->outstream, "%-9s ", get_row_name(lp, i));
535 for(j = 1; j <= lp->columns; j++)
536 fprintf(lp->outstream, "%8g ", get_mat(lp, i, j));
537 if(is_constr_type(lp, i, GE))
538 fprintf(lp->outstream, ">= ");
539 else if(is_constr_type(lp, i, LE))
540 fprintf(lp->outstream, "<= ");
541 else
542 fprintf(lp->outstream, " = ");
543 fprintf(lp->outstream, "%8g", get_rh(lp, i));
544
545 if(is_constr_type(lp, i, GE)) {
546 if(get_rh_upper(lp, i) < lp->infinite)
547 fprintf(lp->outstream, " %s = %8g", "upbo", get_rh_upper(lp, i));
548 }
549 else if(is_constr_type(lp, i, LE)) {
550 if(get_rh_lower(lp, i) > -lp->infinite)
551 fprintf(lp->outstream, " %s = %8g", "lowbo", get_rh_lower(lp, i));
552 }
553 fprintf(lp->outstream, "\n");
554 }
555
556 fprintf(lp->outstream, "Type ");
557 for(i = 1; i <= lp->columns; i++) {
558 if(is_int(lp,i))
559 fprintf(lp->outstream, " Int ");
560 else
561 fprintf(lp->outstream, " Real ");
562 }
563
564 fprintf(lp->outstream, "\nupbo ");
565 for(i = 1; i <= lp->columns; i++)
566 if(get_upbo(lp, i) >= lp->infinite)
567 fprintf(lp->outstream, " Inf ");
568 else
569 fprintf(lp->outstream, "%8g ", get_upbo(lp, i));
570 fprintf(lp->outstream, "\nlowbo ");
571 for(i = 1; i <= lp->columns; i++)
572 if(get_lowbo(lp, i) <= -lp->infinite)
573 fprintf(lp->outstream, " -Inf ");
574 else
575 fprintf(lp->outstream, "%8g ", get_lowbo(lp, i));
576 fprintf(lp->outstream, "\n");
577
578 fflush(lp->outstream);
579 }
580
581 /* Report the scaling factors used; extremely rarely used */
REPORT_scales(lprec * lp)582 void REPORT_scales(lprec *lp)
583 {
584 int i, colMax;
585
586 colMax = lp->columns;
587
588 if(lp->outstream == NULL)
589 return;
590
591 if(lp->scaling_used) {
592 fprintf(lp->outstream, "\nScale factors:\n");
593 for(i = 0; i <= lp->rows + colMax; i++)
594 fprintf(lp->outstream, "%-20s scaled at %g\n",
595 (i <= lp->rows) ? get_row_name(lp, i) : get_col_name(lp, i - lp->rows),
596 (double)lp->scalars[i]);
597 }
598 fflush(lp->outstream);
599 }
600
601 /* Report the traditional tableau corresponding to the current basis */
REPORT_tableau(lprec * lp)602 MYBOOL REPORT_tableau(lprec *lp)
603 {
604 int j, row_nr, *coltarget;
605 REAL *prow = NULL;
606 FILE *stream = lp->outstream;
607
608 if(lp->outstream == NULL)
609 return(FALSE);
610
611 if(!lp->model_is_valid || !has_BFP(lp) ||
612 (get_total_iter(lp) == 0) || (lp->spx_status == NOTRUN)) {
613 lp->spx_status = NOTRUN;
614 return(FALSE);
615 }
616 if(!allocREAL(lp, &prow,lp->sum + 1, TRUE)) {
617 lp->spx_status = NOMEMORY;
618 return(FALSE);
619 }
620
621 fprintf(stream, "\n");
622 fprintf(stream, "Tableau at iter %.0f:\n", (double) get_total_iter(lp));
623
624 for(j = 1; j <= lp->sum; j++)
625 if (!lp->is_basic[j])
626 fprintf(stream, "%15d", (j <= lp->rows ?
627 (j + lp->columns) * ((lp->orig_upbo[j] == 0) ||
628 (is_chsign(lp, j)) ? 1 : -1) : j - lp->rows) *
629 (lp->is_lower[j] ? 1 : -1));
630 fprintf(stream, "\n");
631
632 coltarget = (int *) mempool_obtainVector(lp->workarrays, lp->columns+1, sizeof(*coltarget));
633 if(!get_colIndexA(lp, SCAN_USERVARS+USE_NONBASICVARS, coltarget, FALSE)) {
634 mempool_releaseVector(lp->workarrays, (char *) coltarget, FALSE);
635 return(FALSE);
636 }
637 for(row_nr = 1; (row_nr <= lp->rows + 1); row_nr++) {
638 if (row_nr <= lp->rows)
639 fprintf(stream, "%3d", (lp->var_basic[row_nr] <= lp->rows ?
640 (lp->var_basic[row_nr] + lp->columns) * ((lp->orig_upbo[lp->var_basic [row_nr]] == 0) ||
641 (is_chsign(lp, lp->var_basic[row_nr])) ? 1 : -1) : lp->var_basic[row_nr] - lp->rows) *
642 (lp->is_lower[lp->var_basic [row_nr]] ? 1 : -1));
643 else
644 fprintf(stream, " ");
645 bsolve(lp, row_nr <= lp->rows ? row_nr : 0, prow, NULL, lp->epsmachine*DOUBLEROUND, 1.0);
646 prod_xA(lp, coltarget, prow, NULL, lp->epsmachine, 1.0,
647 prow, NULL, MAT_ROUNDDEFAULT);
648
649 for(j = 1; j <= lp->rows + lp->columns; j++)
650 if (!lp->is_basic[j])
651 fprintf(stream, "%15.7f", prow[j] * (lp->is_lower[j] ? 1 : -1) *
652 (row_nr <= lp->rows ? 1 : -1));
653 fprintf(stream, "%15.7f", lp->rhs[row_nr <= lp->rows ? row_nr : 0] *
654 (double) ((row_nr <= lp->rows) || (is_maxim(lp)) ? 1 : -1));
655 fprintf(stream, "\n");
656 }
657 fflush(stream);
658
659 mempool_releaseVector(lp->workarrays, (char *) coltarget, FALSE);
660 FREE(prow);
661 return(TRUE);
662 }
663
REPORT_constraintinfo(lprec * lp,char * datainfo)664 void REPORT_constraintinfo(lprec *lp, char *datainfo)
665 {
666 int i, tally[ROWCLASS_MAX+1];
667
668 MEMCLEAR(tally, ROWCLASS_MAX+1);
669 for(i = 1; i <= lp->rows; i++)
670 tally[get_constr_class(lp, i)]++;
671
672 if(datainfo != NULL)
673 report(lp, NORMAL, "%s\n", datainfo);
674
675 for(i = 0; i <= ROWCLASS_MAX; i++)
676 if(tally[i] > 0)
677 report(lp, NORMAL, "%-15s %4d\n", get_str_constr_class(lp, i), tally[i]);
678 }
679
REPORT_modelinfo(lprec * lp,MYBOOL doName,char * datainfo)680 void REPORT_modelinfo(lprec *lp, MYBOOL doName, char *datainfo)
681 {
682 if(doName) {
683 report(lp, NORMAL, "\nModel name: '%s' - run #%-5d\n",
684 get_lp_name(lp), lp->solvecount);
685 report(lp, NORMAL, "Objective: %simize(%s)\n",
686 my_if(is_maxim(lp), "Max", "Min"), get_row_name(lp, 0));
687 report(lp, NORMAL, " \n");
688 }
689 if(datainfo != NULL)
690 report(lp, NORMAL, "%s\n", datainfo);
691
692 report(lp, NORMAL, "Model size: %7d constraints, %7d variables, %12d non-zeros.\n",
693 lp->rows, lp->columns, get_nonzeros(lp));
694 if(GUB_count(lp)+SOS_count(lp) > 0)
695 report(lp, NORMAL, "Var-types: %7d integer, %7d semi-cont., %7d SOS.\n",
696 lp->int_vars, lp->sc_vars, lp->sos_vars);
697 report(lp, NORMAL, "Sets: %7d GUB, %7d SOS.\n",
698 GUB_count(lp), SOS_count(lp));
699 }
700
701 /* Save a matrix column subset to a MatrixMarket formatted file,
702 say to export the basis matrix for further numerical analysis.
703 If colndx is NULL, then the full constraint matrix is assumed. */
REPORT_mat_mmsave(lprec * lp,char * filename,int * colndx,MYBOOL includeOF,char * infotext)704 MYBOOL REPORT_mat_mmsave(lprec *lp, char *filename, int *colndx, MYBOOL includeOF, char *infotext)
705 {
706 int n, m, nz, i, j, k, kk;
707 MATrec *mat = lp->matA;
708 MM_typecode matcode;
709 FILE *output; /* = stdout; */
710 MYBOOL ok;
711 REAL *acol = NULL;
712 int *nzlist = NULL;
713
714 /* Open file */
715 ok = (MYBOOL) ((filename == NULL) || ((output = fopen(filename,"w")) != NULL));
716 if(!ok)
717 return(ok);
718 if((filename == NULL) && (lp->outstream != NULL))
719 output = lp->outstream;
720
721 /* Compute column and non-zero counts */
722 if(colndx == lp->var_basic) {
723 if(!lp->basis_valid)
724 return( FALSE );
725 m = lp->rows;
726 }
727 else if(colndx != NULL)
728 m = colndx[0];
729 else
730 m = lp->columns;
731 n = lp->rows;
732 nz = 0;
733
734 for(j = 1; j <= m; j++) {
735 k = (colndx == NULL ? n + j : colndx[j]);
736 if(k > n) {
737 k -= lp->rows;
738 nz += mat_collength(mat, k);
739 if(includeOF && is_OF_nz(lp, k))
740 nz++;
741 }
742 else
743 nz++;
744 }
745 kk = 0;
746 if(includeOF) {
747 n++; /* Row count */
748 kk++; /* Row index offset */
749 }
750
751 /* Initialize */
752 mm_initialize_typecode(&matcode);
753 mm_set_matrix(&matcode);
754 mm_set_coordinate(&matcode);
755 mm_set_real(&matcode);
756
757 mm_write_banner(output, matcode);
758 mm_write_mtx_crd_size(output, n+kk, m, nz+(colndx == lp->var_basic ? 1 : 0));
759
760 /* Allocate working arrays for sparse column storage */
761 allocREAL(lp, &acol, n+2, FALSE);
762 allocINT(lp, &nzlist, n+2, FALSE);
763
764 /* Write the matrix non-zero values column-by-column.
765 NOTE: matrixMarket files use 1-based indeces,
766 i.e. first row of a vector has index 1, not 0. */
767 if(infotext != NULL) {
768 fprintf(output, "%%\n");
769 fprintf(output, "%% %s\n", infotext);
770 fprintf(output, "%%\n");
771 }
772 if(includeOF && (colndx == lp->var_basic))
773 fprintf(output, "%d %d %g\n", 1, 1, 1.0);
774 for(j = 1; j <= m; j++) {
775 k = (colndx == NULL ? lp->rows + j : colndx[j]);
776 if(k == 0)
777 continue;
778 nz = obtain_column(lp, k, acol, nzlist, NULL);
779 for(i = 1; i <= nz; i++) {
780 if(!includeOF && (nzlist[i] == 0))
781 continue;
782 fprintf(output, "%d %d %g\n", nzlist[i]+kk, j+kk, acol[i]);
783 }
784 }
785 fprintf(output, "%% End of MatrixMarket file\n");
786
787 /* Finish */
788 FREE(acol);
789 FREE(nzlist);
790 fclose(output);
791
792 return(ok);
793 }
794
795