1 /* prrngs.c (print sensitivity analysis report) */
2
3 /***********************************************************************
4 * This code is part of GLPK (GNU Linear Programming Kit).
5 * Copyright (C) 2009-2016 Free Software Foundation, Inc.
6 * Written by Andrew Makhorin <mao@gnu.org>.
7 *
8 * GLPK is free software: you can redistribute it and/or modify it
9 * under the terms of the GNU General Public License as published by
10 * the Free Software Foundation, either version 3 of the License, or
11 * (at your option) any later version.
12 *
13 * GLPK is distributed in the hope that it will be useful, but WITHOUT
14 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
16 * License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with GLPK. If not, see <http://www.gnu.org/licenses/>.
20 ***********************************************************************/
21
22 #include "env.h"
23 #include "prob.h"
24
25 #define xfprintf glp_format
26
format(char buf[13+1],double x)27 static char *format(char buf[13+1], double x)
28 { /* format floating-point number in MPS/360-like style */
29 if (x == -DBL_MAX)
30 strcpy(buf, " -Inf");
31 else if (x == +DBL_MAX)
32 strcpy(buf, " +Inf");
33 else if (fabs(x) <= 999999.99998)
34 { sprintf(buf, "%13.5f", x);
35 #if 1
36 if (strcmp(buf, " 0.00000") == 0 ||
37 strcmp(buf, " -0.00000") == 0)
38 strcpy(buf, " . ");
39 else if (memcmp(buf, " 0.", 8) == 0)
40 memcpy(buf, " .", 8);
41 else if (memcmp(buf, " -0.", 8) == 0)
42 memcpy(buf, " -.", 8);
43 #endif
44 }
45 else
46 sprintf(buf, "%13.6g", x);
47 return buf;
48 }
49
glp_print_ranges(glp_prob * P,int len,const int list[],int flags,const char * fname)50 int glp_print_ranges(glp_prob *P, int len, const int list[],
51 int flags, const char *fname)
52 { /* print sensitivity analysis report */
53 glp_file *fp = NULL;
54 GLPROW *row;
55 GLPCOL *col;
56 int m, n, pass, k, t, numb, type, stat, var1, var2, count, page,
57 ret;
58 double lb, ub, slack, coef, prim, dual, value1, value2, coef1,
59 coef2, obj1, obj2;
60 const char *name, *limit;
61 char buf[13+1];
62 /* sanity checks */
63 #if 0 /* 04/IV-2016 */
64 if (P == NULL || P->magic != GLP_PROB_MAGIC)
65 xerror("glp_print_ranges: P = %p; invalid problem object\n",
66 P);
67 #endif
68 m = P->m, n = P->n;
69 if (len < 0)
70 xerror("glp_print_ranges: len = %d; invalid list length\n",
71 len);
72 if (len > 0)
73 { if (list == NULL)
74 xerror("glp_print_ranges: list = %p: invalid parameter\n",
75 list);
76 for (t = 1; t <= len; t++)
77 { k = list[t];
78 if (!(1 <= k && k <= m+n))
79 xerror("glp_print_ranges: list[%d] = %d; row/column numb"
80 "er out of range\n", t, k);
81 }
82 }
83 if (flags != 0)
84 xerror("glp_print_ranges: flags = %d; invalid parameter\n",
85 flags);
86 if (fname == NULL)
87 xerror("glp_print_ranges: fname = %p; invalid parameter\n",
88 fname);
89 if (glp_get_status(P) != GLP_OPT)
90 { xprintf("glp_print_ranges: optimal basic solution required\n");
91 ret = 1;
92 goto done;
93 }
94 if (!glp_bf_exists(P))
95 { xprintf("glp_print_ranges: basis factorization required\n");
96 ret = 2;
97 goto done;
98 }
99 /* start reporting */
100 xprintf("Write sensitivity analysis report to '%s'...\n", fname);
101 fp = glp_open(fname, "w");
102 if (fp == NULL)
103 { xprintf("Unable to create '%s' - %s\n", fname, get_err_msg());
104 ret = 3;
105 goto done;
106 }
107 page = count = 0;
108 for (pass = 1; pass <= 2; pass++)
109 for (t = 1; t <= (len == 0 ? m+n : len); t++)
110 { if (t == 1) count = 0;
111 k = (len == 0 ? t : list[t]);
112 if (pass == 1 && k > m || pass == 2 && k <= m)
113 continue;
114 if (count == 0)
115 { xfprintf(fp, "GLPK %-4s - SENSITIVITY ANALYSIS REPORT%73sPa"
116 "ge%4d\n", glp_version(), "", ++page);
117 xfprintf(fp, "\n");
118 xfprintf(fp, "%-12s%s\n", "Problem:",
119 P->name == NULL ? "" : P->name);
120 xfprintf(fp, "%-12s%s%s%.10g (%s)\n", "Objective:",
121 P->obj == NULL ? "" : P->obj,
122 P->obj == NULL ? "" : " = ", P->obj_val,
123 P->dir == GLP_MIN ? "MINimum" :
124 P->dir == GLP_MAX ? "MAXimum" : "???");
125 xfprintf(fp, "\n");
126 xfprintf(fp, "%6s %-12s %2s %13s %13s %13s %13s %13s %13s "
127 "%s\n", "No.", pass == 1 ? "Row name" : "Column name",
128 "St", "Activity", pass == 1 ? "Slack" : "Obj coef",
129 "Lower bound", "Activity", "Obj coef", "Obj value at",
130 "Limiting");
131 xfprintf(fp, "%6s %-12s %2s %13s %13s %13s %13s %13s %13s "
132 "%s\n", "", "", "", "", "Marginal", "Upper bound",
133 "range", "range", "break point", "variable");
134 xfprintf(fp, "------ ------------ -- ------------- --------"
135 "----- ------------- ------------- ------------- ------"
136 "------- ------------\n");
137 }
138 if (pass == 1)
139 { numb = k;
140 xassert(1 <= numb && numb <= m);
141 row = P->row[numb];
142 name = row->name;
143 type = row->type;
144 lb = glp_get_row_lb(P, numb);
145 ub = glp_get_row_ub(P, numb);
146 coef = 0.0;
147 stat = row->stat;
148 prim = row->prim;
149 if (type == GLP_FR)
150 slack = - prim;
151 else if (type == GLP_LO)
152 slack = lb - prim;
153 else if (type == GLP_UP || type == GLP_DB || type == GLP_FX)
154 slack = ub - prim;
155 dual = row->dual;
156 }
157 else
158 { numb = k - m;
159 xassert(1 <= numb && numb <= n);
160 col = P->col[numb];
161 name = col->name;
162 lb = glp_get_col_lb(P, numb);
163 ub = glp_get_col_ub(P, numb);
164 coef = col->coef;
165 stat = col->stat;
166 prim = col->prim;
167 slack = 0.0;
168 dual = col->dual;
169 }
170 if (stat != GLP_BS)
171 { glp_analyze_bound(P, k, &value1, &var1, &value2, &var2);
172 if (stat == GLP_NF)
173 coef1 = coef2 = coef;
174 else if (stat == GLP_NS)
175 coef1 = -DBL_MAX, coef2 = +DBL_MAX;
176 else if (stat == GLP_NL && P->dir == GLP_MIN ||
177 stat == GLP_NU && P->dir == GLP_MAX)
178 coef1 = coef - dual, coef2 = +DBL_MAX;
179 else
180 coef1 = -DBL_MAX, coef2 = coef - dual;
181 if (value1 == -DBL_MAX)
182 { if (dual < -1e-9)
183 obj1 = +DBL_MAX;
184 else if (dual > +1e-9)
185 obj1 = -DBL_MAX;
186 else
187 obj1 = P->obj_val;
188 }
189 else
190 obj1 = P->obj_val + dual * (value1 - prim);
191 if (value2 == +DBL_MAX)
192 { if (dual < -1e-9)
193 obj2 = -DBL_MAX;
194 else if (dual > +1e-9)
195 obj2 = +DBL_MAX;
196 else
197 obj2 = P->obj_val;
198 }
199 else
200 obj2 = P->obj_val + dual * (value2 - prim);
201 }
202 else
203 { glp_analyze_coef(P, k, &coef1, &var1, &value1, &coef2,
204 &var2, &value2);
205 if (coef1 == -DBL_MAX)
206 { if (prim < -1e-9)
207 obj1 = +DBL_MAX;
208 else if (prim > +1e-9)
209 obj1 = -DBL_MAX;
210 else
211 obj1 = P->obj_val;
212 }
213 else
214 obj1 = P->obj_val + (coef1 - coef) * prim;
215 if (coef2 == +DBL_MAX)
216 { if (prim < -1e-9)
217 obj2 = -DBL_MAX;
218 else if (prim > +1e-9)
219 obj2 = +DBL_MAX;
220 else
221 obj2 = P->obj_val;
222 }
223 else
224 obj2 = P->obj_val + (coef2 - coef) * prim;
225 }
226 /*** first line ***/
227 /* row/column number */
228 xfprintf(fp, "%6d", numb);
229 /* row/column name */
230 xfprintf(fp, " %-12.12s", name == NULL ? "" : name);
231 if (name != NULL && strlen(name) > 12)
232 xfprintf(fp, "%s\n%6s %12s", name+12, "", "");
233 /* row/column status */
234 xfprintf(fp, " %2s",
235 stat == GLP_BS ? "BS" : stat == GLP_NL ? "NL" :
236 stat == GLP_NU ? "NU" : stat == GLP_NF ? "NF" :
237 stat == GLP_NS ? "NS" : "??");
238 /* row/column activity */
239 xfprintf(fp, " %s", format(buf, prim));
240 /* row slack, column objective coefficient */
241 xfprintf(fp, " %s", format(buf, k <= m ? slack : coef));
242 /* row/column lower bound */
243 xfprintf(fp, " %s", format(buf, lb));
244 /* row/column activity range */
245 xfprintf(fp, " %s", format(buf, value1));
246 /* row/column objective coefficient range */
247 xfprintf(fp, " %s", format(buf, coef1));
248 /* objective value at break point */
249 xfprintf(fp, " %s", format(buf, obj1));
250 /* limiting variable name */
251 if (var1 != 0)
252 { if (var1 <= m)
253 limit = glp_get_row_name(P, var1);
254 else
255 limit = glp_get_col_name(P, var1 - m);
256 if (limit != NULL)
257 xfprintf(fp, " %s", limit);
258 }
259 xfprintf(fp, "\n");
260 /*** second line ***/
261 xfprintf(fp, "%6s %-12s %2s %13s", "", "", "", "");
262 /* row/column reduced cost */
263 xfprintf(fp, " %s", format(buf, dual));
264 /* row/column upper bound */
265 xfprintf(fp, " %s", format(buf, ub));
266 /* row/column activity range */
267 xfprintf(fp, " %s", format(buf, value2));
268 /* row/column objective coefficient range */
269 xfprintf(fp, " %s", format(buf, coef2));
270 /* objective value at break point */
271 xfprintf(fp, " %s", format(buf, obj2));
272 /* limiting variable name */
273 if (var2 != 0)
274 { if (var2 <= m)
275 limit = glp_get_row_name(P, var2);
276 else
277 limit = glp_get_col_name(P, var2 - m);
278 if (limit != NULL)
279 xfprintf(fp, " %s", limit);
280 }
281 xfprintf(fp, "\n");
282 xfprintf(fp, "\n");
283 /* print 10 items per page */
284 count = (count + 1) % 10;
285 }
286 xfprintf(fp, "End of report\n");
287 #if 0 /* FIXME */
288 xfflush(fp);
289 #endif
290 if (glp_ioerr(fp))
291 { xprintf("Write error on '%s' - %s\n", fname, get_err_msg());
292 ret = 4;
293 goto done;
294 }
295 ret = 0;
296 done: if (fp != NULL) glp_close(fp);
297 return ret;
298 }
299
300 /* eof */
301