Mercurial > octave
annotate src/DLD-FUNCTIONS/__glpk__.cc @ 8036:854683691d7a
fix invalid memory read in glpk
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Tue, 19 Aug 2008 13:59:30 -0400 |
parents | a1dbe9d80eee |
children | e3e3d12364b0 |
rev | line source |
---|---|
5234 | 1 /* |
2 | |
7017 | 3 Copyright (C) 2005, 2006, 2007 Nicolo' Giorgetti |
5234 | 4 |
5 This file is part of Octave. | |
6 | |
7 Octave is free software; you can redistribute it and/or modify it | |
8 under the terms of the GNU General Public License as published by the | |
7016 | 9 Free Software Foundation; either version 3 of the License, or (at your |
10 option) any later version. | |
5234 | 11 |
12 Octave is distributed in the hope that it will be useful, but WITHOUT | |
13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |
15 for more details. | |
16 | |
17 You should have received a copy of the GNU General Public License | |
7016 | 18 along with Octave; see the file COPYING. If not, see |
19 <http://www.gnu.org/licenses/>. | |
5234 | 20 |
21 */ | |
22 | |
23 #ifdef HAVE_CONFIG_H | |
24 #include <config.h> | |
25 #endif | |
5232 | 26 |
27 #include <cfloat> | |
28 #include <csetjmp> | |
29 #include <ctime> | |
30 | |
5234 | 31 #include "defun-dld.h" |
32 #include "error.h" | |
5235 | 33 #include "gripes.h" |
34 #include "oct-map.h" | |
5234 | 35 #include "oct-obj.h" |
5235 | 36 #include "pager.h" |
37 | |
38 #if defined (HAVE_GLPK) | |
5232 | 39 |
6333 | 40 extern "C" |
41 { | |
6804 | 42 #if defined (HAVE_GLPK_GLPK_H) |
43 #include <glpk/glpk.h> | |
44 #else | |
5234 | 45 #include <glpk.h> |
6804 | 46 #endif |
6333 | 47 |
6381 | 48 #ifdef GLPK_PRE_4_14 |
6333 | 49 |
6287 | 50 #ifndef _GLPLIB_H |
51 #include <glplib.h> | |
52 #endif | |
53 #ifndef lib_set_fault_hook | |
54 #define lib_set_fault_hook lib_fault_hook | |
55 #endif | |
56 #ifndef lib_set_print_hook | |
57 #define lib_set_print_hook lib_print_hook | |
58 #endif | |
5232 | 59 |
6333 | 60 #else |
61 | |
62 void _glp_lib_print_hook (int (*func)(void *info, char *buf), void *info); | |
63 void _glp_lib_fault_hook (int (*func)(void *info, char *buf), void *info); | |
64 | |
65 #endif | |
6472 | 66 } |
6333 | 67 |
5232 | 68 #define NIntP 17 |
69 #define NRealP 10 | |
70 | |
5234 | 71 int lpxIntParam[NIntP] = { |
5240 | 72 0, |
5234 | 73 1, |
74 0, | |
75 1, | |
76 0, | |
77 -1, | |
78 0, | |
79 200, | |
80 1, | |
81 2, | |
82 0, | |
83 1, | |
84 0, | |
85 0, | |
86 2, | |
87 2, | |
88 1 | |
5232 | 89 }; |
90 | |
5234 | 91 int IParam[NIntP] = { |
92 LPX_K_MSGLEV, | |
93 LPX_K_SCALE, | |
94 LPX_K_DUAL, | |
95 LPX_K_PRICE, | |
96 LPX_K_ROUND, | |
97 LPX_K_ITLIM, | |
98 LPX_K_ITCNT, | |
99 LPX_K_OUTFRQ, | |
100 LPX_K_MPSINFO, | |
101 LPX_K_MPSOBJ, | |
102 LPX_K_MPSORIG, | |
103 LPX_K_MPSWIDE, | |
104 LPX_K_MPSFREE, | |
105 LPX_K_MPSSKIP, | |
106 LPX_K_BRANCH, | |
107 LPX_K_BTRACK, | |
108 LPX_K_PRESOL | |
5232 | 109 }; |
110 | |
111 | |
5234 | 112 double lpxRealParam[NRealP] = { |
113 0.07, | |
114 1e-7, | |
115 1e-7, | |
116 1e-9, | |
117 -DBL_MAX, | |
118 DBL_MAX, | |
119 -1.0, | |
120 0.0, | |
121 1e-6, | |
122 1e-7 | |
5232 | 123 }; |
124 | |
5234 | 125 int RParam[NRealP] = { |
126 LPX_K_RELAX, | |
127 LPX_K_TOLBND, | |
128 LPX_K_TOLDJ, | |
129 LPX_K_TOLPIV, | |
130 LPX_K_OBJLL, | |
131 LPX_K_OBJUL, | |
132 LPX_K_TMLIM, | |
133 LPX_K_OUTDLY, | |
134 LPX_K_TOLINT, | |
135 LPX_K_TOLOBJ | |
5232 | 136 }; |
137 | |
5241 | 138 static jmp_buf mark; //-- Address for long jump to jump to |
5232 | 139 |
5234 | 140 int |
5235 | 141 glpk_fault_hook (void * /* info */, char *msg) |
5234 | 142 { |
5240 | 143 error ("CRITICAL ERROR in GLPK: %s", msg); |
5234 | 144 longjmp (mark, -1); |
5232 | 145 } |
146 | |
5234 | 147 int |
5235 | 148 glpk_print_hook (void * /* info */, char *msg) |
5232 | 149 { |
5241 | 150 message (0, "%s", msg); |
5234 | 151 return 1; |
5232 | 152 } |
153 | |
5234 | 154 int |
155 glpk (int sense, int n, int m, double *c, int nz, int *rn, int *cn, | |
156 double *a, double *b, char *ctype, int *freeLB, double *lb, | |
157 int *freeUB, double *ub, int *vartype, int isMIP, int lpsolver, | |
158 int save_pb, double *xmin, double *fmin, double *status, | |
159 double *lambda, double *redcosts, double *time, double *mem) | |
5232 | 160 { |
5240 | 161 int errnum; |
5234 | 162 int typx = 0; |
163 int method; | |
164 | |
165 clock_t t_start = clock(); | |
166 | |
6381 | 167 #ifdef GLPK_PRE_4_14 |
6333 | 168 lib_set_fault_hook (0, glpk_fault_hook); |
169 #else | |
170 _glp_lib_fault_hook (glpk_fault_hook, 0); | |
171 #endif | |
5232 | 172 |
5234 | 173 if (lpxIntParam[0] > 1) |
6381 | 174 #ifdef GLPK_PRE_4_14 |
6333 | 175 lib_set_print_hook (0, glpk_print_hook); |
176 #else | |
177 _glp_lib_print_hook (glpk_print_hook, 0); | |
178 #endif | |
5234 | 179 |
180 LPX *lp = lpx_create_prob (); | |
181 | |
5232 | 182 |
5234 | 183 //-- Set the sense of optimization |
184 if (sense == 1) | |
185 lpx_set_obj_dir (lp, LPX_MIN); | |
186 else | |
187 lpx_set_obj_dir (lp, LPX_MAX); | |
188 | |
189 //-- If the problem has integer structural variables switch to MIP | |
190 if (isMIP) | |
191 lpx_set_class (lp, LPX_MIP); | |
5232 | 192 |
5234 | 193 lpx_add_cols (lp, n); |
194 for (int i = 0; i < n; i++) | |
195 { | |
5232 | 196 //-- Define type of the structural variables |
5234 | 197 if (! freeLB[i] && ! freeUB[i]) |
6898 | 198 { |
199 if (lb[i] != ub[i]) | |
200 lpx_set_col_bnds (lp, i+1, LPX_DB, lb[i], ub[i]); | |
201 else | |
202 lpx_set_col_bnds (lp, i+1, LPX_FX, lb[i], ub[i]); | |
203 } | |
5234 | 204 else |
205 { | |
206 if (! freeLB[i] && freeUB[i]) | |
207 lpx_set_col_bnds (lp, i+1, LPX_LO, lb[i], ub[i]); | |
208 else | |
209 { | |
210 if (freeLB[i] && ! freeUB[i]) | |
211 lpx_set_col_bnds (lp, i+1, LPX_UP, lb[i], ub[i]); | |
212 else | |
213 lpx_set_col_bnds (lp, i+1, LPX_FR, lb[i], ub[i]); | |
214 } | |
215 } | |
216 | |
217 // -- Set the objective coefficient of the corresponding | |
5232 | 218 // -- structural variable. No constant term is assumed. |
5234 | 219 lpx_set_obj_coef(lp,i+1,c[i]); |
5232 | 220 |
5234 | 221 if (isMIP) |
222 lpx_set_col_kind (lp, i+1, vartype[i]); | |
223 } | |
224 | |
225 lpx_add_rows (lp, m); | |
226 | |
227 for (int i = 0; i < m; i++) | |
228 { | |
5232 | 229 /* If the i-th row has no lower bound (types F,U), the |
230 corrispondent parameter will be ignored. | |
231 If the i-th row has no upper bound (types F,L), the corrispondent | |
232 parameter will be ignored. | |
233 If the i-th row is of S type, the i-th LB is used, but | |
234 the i-th UB is ignored. | |
235 */ | |
5234 | 236 |
237 switch (ctype[i]) | |
238 { | |
239 case 'F': | |
240 typx = LPX_FR; | |
241 break; | |
242 | |
243 case 'U': | |
244 typx = LPX_UP; | |
245 break; | |
246 | |
247 case 'L': | |
248 typx = LPX_LO; | |
249 break; | |
5232 | 250 |
5234 | 251 case 'S': |
252 typx = LPX_FX; | |
253 break; | |
5232 | 254 |
5234 | 255 case 'D': |
256 typx = LPX_DB; | |
257 break; | |
258 } | |
259 | |
260 lpx_set_row_bnds (lp, i+1, typx, b[i], b[i]); | |
261 | |
262 } | |
263 | |
264 lpx_load_matrix (lp, nz, rn, cn, a); | |
5232 | 265 |
5234 | 266 if (save_pb) |
267 { | |
6484 | 268 static char tmp[] = "outpb.lp"; |
269 if (lpx_write_cpxlp (lp, tmp) != 0) | |
5234 | 270 { |
5240 | 271 error ("__glpk__: unable to write problem"); |
5234 | 272 longjmp (mark, -1); |
273 } | |
274 } | |
275 | |
276 //-- scale the problem data (if required) | |
277 //-- if (scale && (!presol || method == 1)) lpx_scale_prob(lp); | |
278 //-- LPX_K_SCALE=IParam[1] LPX_K_PRESOL=IParam[16] | |
279 if (lpxIntParam[1] && (! lpxIntParam[16] || lpsolver != 1)) | |
280 lpx_scale_prob (lp); | |
281 | |
282 //-- build advanced initial basis (if required) | |
283 if (lpsolver == 1 && ! lpxIntParam[16]) | |
284 lpx_adv_basis (lp); | |
285 | |
286 for(int i = 0; i < NIntP; i++) | |
287 lpx_set_int_parm (lp, IParam[i], lpxIntParam[i]); | |
5232 | 288 |
5234 | 289 for (int i = 0; i < NRealP; i++) |
290 lpx_set_real_parm (lp, RParam[i], lpxRealParam[i]); | |
291 | |
292 if (lpsolver == 1) | |
293 method = 'S'; | |
294 else | |
295 method = 'T'; | |
5232 | 296 |
5234 | 297 switch (method) |
298 { | |
299 case 'S': | |
300 { | |
301 if (isMIP) | |
302 { | |
303 method = 'I'; | |
5240 | 304 errnum = lpx_simplex (lp); |
305 errnum = lpx_integer (lp); | |
5234 | 306 } |
307 else | |
5240 | 308 errnum = lpx_simplex(lp); |
5234 | 309 } |
5232 | 310 break; |
5234 | 311 |
312 case 'T': | |
5240 | 313 errnum = lpx_interior(lp); |
5234 | 314 break; |
315 | |
316 default: | |
6381 | 317 #ifdef GLPK_PRE_4_14 |
5234 | 318 insist (method != method); |
6333 | 319 #else |
6484 | 320 static char tmp[] = "method != method"; |
321 glpk_fault_hook (0, tmp); | |
6333 | 322 #endif |
5234 | 323 } |
324 | |
5240 | 325 /* errnum assumes the following results: |
326 errnum = 0 <=> No errors | |
327 errnum = 1 <=> Iteration limit exceeded. | |
328 errnum = 2 <=> Numerical problems with basis matrix. | |
5234 | 329 */ |
5240 | 330 if (errnum == LPX_E_OK) |
5234 | 331 { |
332 if (isMIP) | |
333 { | |
5241 | 334 *status = lpx_mip_status (lp); |
5234 | 335 *fmin = lpx_mip_obj_val (lp); |
336 } | |
337 else | |
338 { | |
339 if (lpsolver == 1) | |
340 { | |
5241 | 341 *status = lpx_get_status (lp); |
5234 | 342 *fmin = lpx_get_obj_val (lp); |
343 } | |
344 else | |
345 { | |
5241 | 346 *status = lpx_ipt_status (lp); |
5234 | 347 *fmin = lpx_ipt_obj_val (lp); |
348 } | |
349 } | |
5232 | 350 |
5234 | 351 if (isMIP) |
352 { | |
353 for (int i = 0; i < n; i++) | |
354 xmin[i] = lpx_mip_col_val (lp, i+1); | |
355 } | |
356 else | |
357 { | |
358 /* Primal values */ | |
359 for (int i = 0; i < n; i++) | |
360 { | |
361 if (lpsolver == 1) | |
362 xmin[i] = lpx_get_col_prim (lp, i+1); | |
363 else | |
364 xmin[i] = lpx_ipt_col_prim (lp, i+1); | |
365 } | |
5232 | 366 |
5234 | 367 /* Dual values */ |
368 for (int i = 0; i < m; i++) | |
369 { | |
370 if (lpsolver == 1) | |
371 lambda[i] = lpx_get_row_dual (lp, i+1); | |
372 else | |
373 lambda[i] = lpx_ipt_row_dual (lp, i+1); | |
374 } | |
375 | |
376 /* Reduced costs */ | |
377 for (int i = 0; i < lpx_get_num_cols (lp); i++) | |
378 { | |
379 if (lpsolver == 1) | |
380 redcosts[i] = lpx_get_col_dual (lp, i+1); | |
381 else | |
382 redcosts[i] = lpx_ipt_col_dual (lp, i+1); | |
383 } | |
384 } | |
385 | |
5241 | 386 *time = (clock () - t_start) / CLOCKS_PER_SEC; |
6333 | 387 |
6381 | 388 #ifdef GLPK_PRE_4_14 |
5241 | 389 *mem = (lib_env_ptr () -> mem_tpeak); |
6333 | 390 #else |
391 *mem = 0; | |
392 #endif | |
5234 | 393 |
394 lpx_delete_prob (lp); | |
395 return 0; | |
396 } | |
397 | |
398 lpx_delete_prob (lp); | |
399 | |
5241 | 400 *status = errnum; |
5234 | 401 |
5240 | 402 return errnum; |
5232 | 403 } |
404 | |
5235 | 405 #endif |
406 | |
5240 | 407 #define OCTAVE_GLPK_GET_REAL_PARAM(NAME, IDX) \ |
408 do \ | |
409 { \ | |
410 if (PARAM.contains (NAME)) \ | |
411 { \ | |
412 Cell tmp = PARAM.contents (NAME); \ | |
413 \ | |
414 if (! tmp.is_empty ()) \ | |
415 { \ | |
416 lpxRealParam[IDX] = tmp(0).scalar_value (); \ | |
417 \ | |
418 if (error_state) \ | |
419 { \ | |
420 error ("glpk: invalid value in param." NAME); \ | |
421 return retval; \ | |
422 } \ | |
423 } \ | |
424 else \ | |
425 { \ | |
426 error ("glpk: invalid value in param." NAME); \ | |
427 return retval; \ | |
428 } \ | |
429 } \ | |
430 } \ | |
431 while (0) | |
432 | |
433 #define OCTAVE_GLPK_GET_INT_PARAM(NAME, VAL) \ | |
434 do \ | |
435 { \ | |
436 if (PARAM.contains (NAME)) \ | |
437 { \ | |
438 Cell tmp = PARAM.contents (NAME); \ | |
439 \ | |
440 if (! tmp.is_empty ()) \ | |
441 { \ | |
442 VAL = tmp(0).int_value (); \ | |
443 \ | |
444 if (error_state) \ | |
445 { \ | |
446 error ("glpk: invalid value in param." NAME); \ | |
447 return retval; \ | |
448 } \ | |
449 } \ | |
450 else \ | |
451 { \ | |
452 error ("glpk: invalid value in param." NAME); \ | |
453 return retval; \ | |
454 } \ | |
455 } \ | |
456 } \ | |
457 while (0) | |
458 | |
5235 | 459 DEFUN_DLD (__glpk__, args, , |
5245 | 460 "-*- texinfo -*-\n\ |
461 @deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\n\ | |
6945 | 462 Undocumented internal function.\n\ |
5245 | 463 @end deftypefn") |
5232 | 464 { |
5234 | 465 // The list of values to return. See the declaration in oct-obj.h |
466 octave_value_list retval; | |
5232 | 467 |
5235 | 468 #if defined (HAVE_GLPK) |
469 | |
5234 | 470 int nrhs = args.length (); |
471 | |
5240 | 472 if (nrhs != 9) |
5234 | 473 { |
5823 | 474 print_usage (); |
5232 | 475 return retval; |
5234 | 476 } |
477 | |
5237 | 478 //-- 1nd Input. A column array containing the objective function |
479 //-- coefficients. | |
5880 | 480 volatile int mrowsc = args(0).rows(); |
5234 | 481 |
5237 | 482 Matrix C (args(0).matrix_value ()); |
5240 | 483 |
484 if (error_state) | |
485 { | |
486 error ("__glpk__: invalid value of C"); | |
487 return retval; | |
488 } | |
489 | |
5234 | 490 double *c = C.fortran_vec (); |
5455 | 491 Array<int> rn; |
492 Array<int> cn; | |
493 ColumnVector a; | |
5515 | 494 volatile int mrowsA; |
5455 | 495 volatile int nz = 0; |
5232 | 496 |
5237 | 497 //-- 2nd Input. A matrix containing the constraints coefficients. |
5234 | 498 // If matrix A is NOT a sparse matrix |
5631 | 499 if (args(1).is_sparse_type ()) |
5455 | 500 { |
5603 | 501 SparseMatrix A = args(1).sparse_matrix_value (); // get the sparse matrix |
5455 | 502 |
503 if (error_state) | |
504 { | |
505 error ("__glpk__: invalid value of A"); | |
506 return retval; | |
507 } | |
5234 | 508 |
5455 | 509 mrowsA = A.rows (); |
510 octave_idx_type Anc = A.cols (); | |
5604 | 511 octave_idx_type Anz = A.nzmax (); |
5455 | 512 rn.resize (Anz+1); |
513 cn.resize (Anz+1); | |
514 a.resize (Anz+1, 0.0); | |
515 | |
516 if (Anc != mrowsc) | |
517 { | |
518 error ("__glpk__: invalid value of A"); | |
519 return retval; | |
520 } | |
521 | |
522 for (octave_idx_type j = 0; j < Anc; j++) | |
523 for (octave_idx_type i = A.cidx(j); i < A.cidx(j+1); i++) | |
524 { | |
525 nz++; | |
526 rn(nz) = A.ridx(i) + 1; | |
527 cn(nz) = j + 1; | |
5603 | 528 a(nz) = A.data(i); |
5455 | 529 } |
530 } | |
5631 | 531 else |
532 { | |
533 Matrix A (args(1).matrix_value ()); // get the matrix | |
534 | |
535 if (error_state) | |
536 { | |
537 error ("__glpk__: invalid value of A"); | |
538 return retval; | |
539 } | |
540 | |
541 mrowsA = A.rows (); | |
542 rn.resize (mrowsA*mrowsc+1); | |
543 cn.resize (mrowsA*mrowsc+1); | |
544 a.resize (mrowsA*mrowsc+1, 0.0); | |
545 | |
546 for (int i = 0; i < mrowsA; i++) | |
547 { | |
548 for (int j = 0; j < mrowsc; j++) | |
549 { | |
550 if (A(i,j) != 0) | |
551 { | |
552 nz++; | |
553 rn(nz) = i + 1; | |
554 cn(nz) = j + 1; | |
555 a(nz) = A(i,j); | |
556 } | |
557 } | |
558 } | |
559 | |
560 } | |
5232 | 561 |
5237 | 562 //-- 3rd Input. A column array containing the right-hand side value |
5234 | 563 // for each constraint in the constraint matrix. |
5237 | 564 Matrix B (args(2).matrix_value ()); |
5240 | 565 |
566 if (error_state) | |
567 { | |
568 error ("__glpk__: invalid value of b"); | |
569 return retval; | |
570 } | |
571 | |
5234 | 572 double *b = B.fortran_vec (); |
573 | |
5237 | 574 //-- 4th Input. An array of length mrowsc containing the lower |
5234 | 575 //-- bound on each of the variables. |
5237 | 576 Matrix LB (args(3).matrix_value ()); |
5240 | 577 |
8036
854683691d7a
fix invalid memory read in glpk
Jaroslav Hajek <highegg@gmail.com>
parents:
7017
diff
changeset
|
578 if (error_state || LB.length () < mrowsc) |
5240 | 579 { |
580 error ("__glpk__: invalid value of lb"); | |
581 return retval; | |
582 } | |
583 | |
5234 | 584 double *lb = LB.fortran_vec (); |
585 | |
586 //-- LB argument, default: Free | |
587 Array<int> freeLB (mrowsc); | |
588 for (int i = 0; i < mrowsc; i++) | |
589 { | |
5455 | 590 if (xisinf (lb[i])) |
5234 | 591 { |
592 freeLB(i) = 1; | |
593 lb[i] = -octave_Inf; | |
594 } | |
595 else | |
596 freeLB(i) = 0; | |
597 } | |
598 | |
5237 | 599 //-- 5th Input. An array of at least length numcols containing the upper |
5234 | 600 //-- bound on each of the variables. |
5237 | 601 Matrix UB (args(4).matrix_value ()); |
5234 | 602 |
8036
854683691d7a
fix invalid memory read in glpk
Jaroslav Hajek <highegg@gmail.com>
parents:
7017
diff
changeset
|
603 if (error_state || UB.length () < mrowsc) |
5240 | 604 { |
605 error ("__glpk__: invalid value of ub"); | |
606 return retval; | |
607 } | |
608 | |
5234 | 609 double *ub = UB.fortran_vec (); |
5232 | 610 |
5234 | 611 Array<int> freeUB (mrowsc); |
612 for (int i = 0; i < mrowsc; i++) | |
613 { | |
5455 | 614 if (xisinf (ub[i])) |
5234 | 615 { |
616 freeUB(i) = 1; | |
617 ub[i] = octave_Inf; | |
618 } | |
619 else | |
620 freeUB(i) = 0; | |
621 } | |
622 | |
5237 | 623 //-- 6th Input. A column array containing the sense of each constraint |
624 //-- in the constraint matrix. | |
625 charMatrix CTYPE (args(5).char_matrix_value ()); | |
5240 | 626 |
627 if (error_state) | |
628 { | |
629 error ("__glpk__: invalid value of ctype"); | |
630 return retval; | |
631 } | |
632 | |
5237 | 633 char *ctype = CTYPE.fortran_vec (); |
634 | |
635 //-- 7th Input. A column array containing the types of the variables. | |
636 charMatrix VTYPE (args(6).char_matrix_value ()); | |
5232 | 637 |
5240 | 638 if (error_state) |
639 { | |
640 error ("__glpk__: invalid value of vtype"); | |
641 return retval; | |
642 } | |
643 | |
5234 | 644 Array<int> vartype (mrowsc); |
5235 | 645 volatile int isMIP = 0; |
5234 | 646 for (int i = 0; i < mrowsc ; i++) |
647 { | |
648 if (VTYPE(i,0) == 'I') | |
649 { | |
650 isMIP = 1; | |
651 vartype(i) = LPX_IV; | |
652 } | |
653 else | |
5235 | 654 vartype(i) = LPX_CV; |
5234 | 655 } |
656 | |
5237 | 657 //-- 8th Input. Sense of optimization. |
658 volatile int sense; | |
659 double SENSE = args(7).scalar_value (); | |
5240 | 660 |
661 if (error_state) | |
662 { | |
663 error ("__glpk__: invalid value of sense"); | |
664 return retval; | |
665 } | |
666 | |
5237 | 667 if (SENSE >= 0) |
668 sense = 1; | |
669 else | |
670 sense = -1; | |
671 | |
5234 | 672 //-- 9th Input. A structure containing the control parameters. |
673 Octave_map PARAM = args(8).map_value (); | |
674 | |
5240 | 675 if (error_state) |
676 { | |
677 error ("__glpk__: invalid value of param"); | |
678 return retval; | |
679 } | |
680 | |
5234 | 681 //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ |
682 //-- Integer parameters | |
683 //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
684 | |
685 //-- Level of messages output by the solver | |
5240 | 686 OCTAVE_GLPK_GET_INT_PARAM ("msglev", lpxIntParam[0]); |
687 if (lpxIntParam[0] < 0 || lpxIntParam[0] > 3) | |
5234 | 688 { |
5240 | 689 error ("__glpk__: param.msglev must be 0 (no output [default]) or 1 (error messages only) or 2 (normal output) or 3 (full output)"); |
690 return retval; | |
691 } | |
5232 | 692 |
5234 | 693 //-- scaling option |
5240 | 694 OCTAVE_GLPK_GET_INT_PARAM ("scale", lpxIntParam[1]); |
695 if (lpxIntParam[1] < 0 || lpxIntParam[1] > 2) | |
5234 | 696 { |
5240 | 697 error ("__glpk__: param.scale must be 0 (no scaling) or 1 (equilibration scaling [default]) or 2 (geometric mean scaling)"); |
698 return retval; | |
5234 | 699 } |
5232 | 700 |
5234 | 701 //-- Dual dimplex option |
5240 | 702 OCTAVE_GLPK_GET_INT_PARAM ("dual", lpxIntParam[2]); |
703 if (lpxIntParam[2] < 0 || lpxIntParam[2] > 1) | |
5234 | 704 { |
5240 | 705 error ("__glpk__: param.dual must be 0 (do NOT use dual simplex [default]) or 1 (use dual simplex)"); |
706 return retval; | |
5234 | 707 } |
5232 | 708 |
5234 | 709 //-- Pricing option |
5240 | 710 OCTAVE_GLPK_GET_INT_PARAM ("price", lpxIntParam[3]); |
711 if (lpxIntParam[3] < 0 || lpxIntParam[3] > 1) | |
5234 | 712 { |
5240 | 713 error ("__glpk__: param.price must be 0 (textbook pricing) or 1 (steepest edge pricing [default])"); |
714 return retval; | |
715 } | |
5232 | 716 |
5234 | 717 //-- Solution rounding option |
5240 | 718 OCTAVE_GLPK_GET_INT_PARAM ("round", lpxIntParam[4]); |
719 if (lpxIntParam[4] < 0 || lpxIntParam[4] > 1) | |
5234 | 720 { |
5240 | 721 error ("__glpk__: param.round must be 0 (report all primal and dual values [default]) or 1 (replace tiny primal and dual values by exact zero)"); |
722 return retval; | |
5234 | 723 } |
724 | |
725 //-- Simplex iterations limit | |
5240 | 726 OCTAVE_GLPK_GET_INT_PARAM ("itlim", lpxIntParam[5]); |
5232 | 727 |
5234 | 728 //-- Simplex iterations count |
5240 | 729 OCTAVE_GLPK_GET_INT_PARAM ("itcnt", lpxIntParam[6]); |
5234 | 730 |
731 //-- Output frequency, in iterations | |
5240 | 732 OCTAVE_GLPK_GET_INT_PARAM ("outfrq", lpxIntParam[7]); |
5234 | 733 |
734 //-- Branching heuristic option | |
5240 | 735 OCTAVE_GLPK_GET_INT_PARAM ("branch", lpxIntParam[14]); |
736 if (lpxIntParam[14] < 0 || lpxIntParam[14] > 2) | |
5234 | 737 { |
5240 | 738 error ("__glpk__: param.branch must be (MIP only) 0 (branch on first variable) or 1 (branch on last variable) or 2 (branch using a heuristic by Driebeck and Tomlin [default]"); |
739 return retval; | |
740 } | |
5232 | 741 |
5234 | 742 //-- Backtracking heuristic option |
5240 | 743 OCTAVE_GLPK_GET_INT_PARAM ("btrack", lpxIntParam[15]); |
744 if (lpxIntParam[15] < 0 || lpxIntParam[15] > 2) | |
5234 | 745 { |
5240 | 746 error ("__glpk__: param.btrack must be (MIP only) 0 (depth first search) or 1 (breadth first search) or 2 (backtrack using the best projection heuristic [default]"); |
747 return retval; | |
748 } | |
5232 | 749 |
5234 | 750 //-- Presolver option |
5240 | 751 OCTAVE_GLPK_GET_INT_PARAM ("presol", lpxIntParam[16]); |
752 if (lpxIntParam[16] < 0 || lpxIntParam[16] > 1) | |
5234 | 753 { |
5240 | 754 error ("__glpk__: param.presol must be 0 (do NOT use LP presolver) or 1 (use LP presolver [default])"); |
755 return retval; | |
5234 | 756 } |
757 | |
5237 | 758 //-- LPsolver option |
759 volatile int lpsolver = 1; | |
5240 | 760 OCTAVE_GLPK_GET_INT_PARAM ("lpsolver", lpsolver); |
761 if (lpsolver < 1 || lpsolver > 2) | |
5237 | 762 { |
5240 | 763 error ("__glpk__: param.lpsolver must be 1 (simplex method) or 2 (interior point method)"); |
764 return retval; | |
5237 | 765 } |
766 | |
767 //-- Save option | |
768 volatile int save_pb = 0; | |
5240 | 769 OCTAVE_GLPK_GET_INT_PARAM ("save", save_pb); |
770 save_pb = save_pb != 0; | |
5237 | 771 |
5234 | 772 //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ |
773 //-- Real parameters | |
774 //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |
775 | |
776 //-- Ratio test option | |
5240 | 777 OCTAVE_GLPK_GET_REAL_PARAM ("relax", 0); |
5232 | 778 |
5234 | 779 //-- Relative tolerance used to check if the current basic solution |
780 //-- is primal feasible | |
5240 | 781 OCTAVE_GLPK_GET_REAL_PARAM ("tolbnd", 1); |
5234 | 782 |
783 //-- Absolute tolerance used to check if the current basic solution | |
784 //-- is dual feasible | |
5240 | 785 OCTAVE_GLPK_GET_REAL_PARAM ("toldj", 2); |
5232 | 786 |
5234 | 787 //-- Relative tolerance used to choose eligible pivotal elements of |
788 //-- the simplex table in the ratio test | |
5240 | 789 OCTAVE_GLPK_GET_REAL_PARAM ("tolpiv", 3); |
5234 | 790 |
5240 | 791 OCTAVE_GLPK_GET_REAL_PARAM ("objll", 4); |
5234 | 792 |
5240 | 793 OCTAVE_GLPK_GET_REAL_PARAM ("objul", 5); |
5232 | 794 |
5240 | 795 OCTAVE_GLPK_GET_REAL_PARAM ("tmlim", 6); |
5234 | 796 |
5240 | 797 OCTAVE_GLPK_GET_REAL_PARAM ("outdly", 7); |
5234 | 798 |
5240 | 799 OCTAVE_GLPK_GET_REAL_PARAM ("tolint", 8); |
5234 | 800 |
5240 | 801 OCTAVE_GLPK_GET_REAL_PARAM ("tolobj", 9); |
5234 | 802 |
803 //-- Assign pointers to the output parameters | |
804 ColumnVector xmin (mrowsc); | |
805 ColumnVector fmin (1); | |
806 ColumnVector status (1); | |
807 ColumnVector lambda (mrowsA); | |
808 ColumnVector redcosts (mrowsc); | |
809 ColumnVector time (1); | |
810 ColumnVector mem (1); | |
5232 | 811 |
5234 | 812 int jmpret = setjmp (mark); |
5235 | 813 |
5234 | 814 if (jmpret == 0) |
5235 | 815 glpk (sense, mrowsc, mrowsA, c, nz, rn.fortran_vec (), |
816 cn.fortran_vec (), a.fortran_vec (), b, ctype, | |
817 freeLB.fortran_vec (), lb, freeUB.fortran_vec (), | |
818 ub, vartype.fortran_vec (), isMIP, lpsolver, | |
819 save_pb, xmin.fortran_vec (), fmin.fortran_vec (), | |
820 status.fortran_vec (), lambda.fortran_vec (), | |
821 redcosts.fortran_vec (), time.fortran_vec (), | |
822 mem.fortran_vec ()); | |
5232 | 823 |
5234 | 824 Octave_map extra; |
825 | |
5238 | 826 if (! isMIP) |
827 { | |
828 extra.assign ("lambda", octave_value (lambda)); | |
829 extra.assign ("redcosts", octave_value (redcosts)); | |
830 } | |
831 | |
5234 | 832 extra.assign ("time", octave_value (time)); |
833 extra.assign ("mem", octave_value (mem)); | |
5232 | 834 |
5234 | 835 retval(3) = extra; |
836 retval(2) = octave_value(status); | |
837 retval(1) = octave_value(fmin); | |
838 retval(0) = octave_value(xmin); | |
839 | |
5235 | 840 #else |
841 | |
842 gripe_not_supported ("glpk"); | |
843 | |
844 #endif | |
845 | |
5234 | 846 return retval; |
5232 | 847 } |