File: | libinterp/dldfcn/__glpk__.cc |
Location: | line 739, column 25 |
Description: | Function call argument is an uninitialized value |
1 | /* | |||
2 | ||||
3 | Copyright (C) 2005-2013 Nicolo' Giorgetti | |||
4 | Copyright (C) 2013 Sébastien Villemot <sebastien@debian.org> | |||
5 | ||||
6 | This file is part of Octave. | |||
7 | ||||
8 | Octave is free software; you can redistribute it and/or modify it | |||
9 | under the terms of the GNU General Public License as published by the | |||
10 | Free Software Foundation; either version 3 of the License, or (at your | |||
11 | option) any later version. | |||
12 | ||||
13 | Octave is distributed in the hope that it will be useful, but WITHOUT | |||
14 | ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |||
15 | FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |||
16 | for more details. | |||
17 | ||||
18 | You should have received a copy of the GNU General Public License | |||
19 | along with Octave; see the file COPYING. If not, see | |||
20 | <http://www.gnu.org/licenses/>. | |||
21 | ||||
22 | */ | |||
23 | ||||
24 | #ifdef HAVE_CONFIG_H1 | |||
25 | #include <config.h> | |||
26 | #endif | |||
27 | ||||
28 | #include <cfloat> | |||
29 | #include <csetjmp> | |||
30 | #include <ctime> | |||
31 | ||||
32 | #include "lo-ieee.h" | |||
33 | ||||
34 | #include "defun-dld.h" | |||
35 | #include "error.h" | |||
36 | #include "gripes.h" | |||
37 | #include "oct-map.h" | |||
38 | #include "oct-obj.h" | |||
39 | #include "pager.h" | |||
40 | ||||
41 | #if defined (HAVE_GLPK1) | |||
42 | ||||
43 | extern "C" | |||
44 | { | |||
45 | #if defined (HAVE_GLPK_GLPK_H) | |||
46 | #include <glpk/glpk.h> | |||
47 | #else | |||
48 | #include <glpk.h> | |||
49 | #endif | |||
50 | } | |||
51 | ||||
52 | struct control_params | |||
53 | { | |||
54 | int msglev; | |||
55 | int dual; | |||
56 | int price; | |||
57 | int itlim; | |||
58 | int outfrq; | |||
59 | int branch; | |||
60 | int btrack; | |||
61 | int presol; | |||
62 | int rtest; | |||
63 | int tmlim; | |||
64 | int outdly; | |||
65 | double tolbnd; | |||
66 | double toldj; | |||
67 | double tolpiv; | |||
68 | double objll; | |||
69 | double objul; | |||
70 | double tolint; | |||
71 | double tolobj; | |||
72 | }; | |||
73 | ||||
74 | static jmp_buf mark; //-- Address for long jump to jump to | |||
75 | ||||
76 | int | |||
77 | glpk (int sense, int n, int m, double *c, int nz, int *rn, int *cn, | |||
78 | double *a, double *b, char *ctype, int *freeLB, double *lb, | |||
79 | int *freeUB, double *ub, int *vartype, int isMIP, int lpsolver, | |||
80 | int save_pb, int scale, const control_params *par, | |||
81 | double *xmin, double *fmin, int *status, | |||
82 | double *lambda, double *redcosts, double *time) | |||
83 | { | |||
84 | int typx = 0; | |||
85 | int errnum = 0; | |||
86 | ||||
87 | clock_t t_start = clock (); | |||
88 | ||||
89 | glp_prob *lp = glp_create_prob (); | |||
90 | ||||
91 | //-- Set the sense of optimization | |||
92 | if (sense == 1) | |||
93 | glp_set_obj_dir (lp, GLP_MIN1); | |||
94 | else | |||
95 | glp_set_obj_dir (lp, GLP_MAX2); | |||
96 | ||||
97 | glp_add_cols (lp, n); | |||
98 | for (int i = 0; i < n; i++) | |||
99 | { | |||
100 | //-- Define type of the structural variables | |||
101 | if (! freeLB[i] && ! freeUB[i]) | |||
102 | { | |||
103 | if (lb[i] != ub[i]) | |||
104 | glp_set_col_bnds (lp, i+1, GLP_DB4, lb[i], ub[i]); | |||
105 | else | |||
106 | glp_set_col_bnds (lp, i+1, GLP_FX5, lb[i], ub[i]); | |||
107 | } | |||
108 | else | |||
109 | { | |||
110 | if (! freeLB[i] && freeUB[i]) | |||
111 | glp_set_col_bnds (lp, i+1, GLP_LO2, lb[i], ub[i]); | |||
112 | else | |||
113 | { | |||
114 | if (freeLB[i] && ! freeUB[i]) | |||
115 | glp_set_col_bnds (lp, i+1, GLP_UP3, lb[i], ub[i]); | |||
116 | else | |||
117 | glp_set_col_bnds (lp, i+1, GLP_FR1, lb[i], ub[i]); | |||
118 | } | |||
119 | } | |||
120 | ||||
121 | // -- Set the objective coefficient of the corresponding | |||
122 | // -- structural variable. No constant term is assumed. | |||
123 | glp_set_obj_coef(lp,i+1,c[i]); | |||
124 | ||||
125 | if (isMIP) | |||
126 | glp_set_col_kind (lp, i+1, vartype[i]); | |||
127 | } | |||
128 | ||||
129 | glp_add_rows (lp, m); | |||
130 | ||||
131 | for (int i = 0; i < m; i++) | |||
132 | { | |||
133 | /* If the i-th row has no lower bound (types F,U), the | |||
134 | corrispondent parameter will be ignored. | |||
135 | If the i-th row has no upper bound (types F,L), the corrispondent | |||
136 | parameter will be ignored. | |||
137 | If the i-th row is of S type, the i-th LB is used, but | |||
138 | the i-th UB is ignored. | |||
139 | */ | |||
140 | ||||
141 | switch (ctype[i]) | |||
142 | { | |||
143 | case 'F': | |||
144 | typx = GLP_FR1; | |||
145 | break; | |||
146 | ||||
147 | case 'U': | |||
148 | typx = GLP_UP3; | |||
149 | break; | |||
150 | ||||
151 | case 'L': | |||
152 | typx = GLP_LO2; | |||
153 | break; | |||
154 | ||||
155 | case 'S': | |||
156 | typx = GLP_FX5; | |||
157 | break; | |||
158 | ||||
159 | case 'D': | |||
160 | typx = GLP_DB4; | |||
161 | break; | |||
162 | } | |||
163 | ||||
164 | glp_set_row_bnds (lp, i+1, typx, b[i], b[i]); | |||
165 | ||||
166 | } | |||
167 | ||||
168 | glp_load_matrix (lp, nz, rn, cn, a); | |||
169 | ||||
170 | if (save_pb) | |||
171 | { | |||
172 | static char tmp[] = "outpb.lp"; | |||
173 | if (glp_write_lp (lp, NULL__null, tmp) != 0) | |||
174 | { | |||
175 | error ("__glpk__: unable to write problem"); | |||
176 | longjmp (mark, -1); | |||
177 | } | |||
178 | } | |||
179 | ||||
180 | //-- scale the problem data | |||
181 | if (!par->presol || lpsolver != 1) | |||
182 | glp_scale_prob (lp, scale); | |||
183 | ||||
184 | //-- build advanced initial basis (if required) | |||
185 | if (lpsolver == 1 && !par->presol) | |||
186 | glp_adv_basis (lp, 0); | |||
187 | ||||
188 | /* For MIP problems without a presolver, a first pass with glp_simplex | |||
189 | is required */ | |||
190 | if ((!isMIP && lpsolver == 1) | |||
191 | || (isMIP && !par->presol)) | |||
192 | { | |||
193 | glp_smcp smcp; | |||
194 | glp_init_smcp (&smcp); | |||
195 | smcp.msg_lev = par->msglev; | |||
196 | smcp.meth = par->dual; | |||
197 | smcp.pricing = par->price; | |||
198 | smcp.r_test = par->rtest; | |||
199 | smcp.tol_bnd = par->tolbnd; | |||
200 | smcp.tol_dj = par->toldj; | |||
201 | smcp.tol_piv = par->tolpiv; | |||
202 | smcp.obj_ll = par->objll; | |||
203 | smcp.obj_ul = par->objul; | |||
204 | smcp.it_lim = par->itlim; | |||
205 | smcp.tm_lim = par->tmlim; | |||
206 | smcp.out_frq = par->outfrq; | |||
207 | smcp.out_dly = par->outdly; | |||
208 | smcp.presolve = par->presol; | |||
209 | errnum = glp_simplex (lp, &smcp); | |||
210 | } | |||
211 | ||||
212 | if (isMIP) | |||
213 | { | |||
214 | glp_iocp iocp; | |||
215 | glp_init_iocp (&iocp); | |||
216 | iocp.msg_lev = par->msglev; | |||
217 | iocp.br_tech = par->branch; | |||
218 | iocp.bt_tech = par->btrack; | |||
219 | iocp.tol_int = par->tolint; | |||
220 | iocp.tol_obj = par->tolobj; | |||
221 | iocp.tm_lim = par->tmlim; | |||
222 | iocp.out_frq = par->outfrq; | |||
223 | iocp.out_dly = par->outdly; | |||
224 | iocp.presolve = par->presol; | |||
225 | errnum = glp_intopt (lp, &iocp); | |||
226 | } | |||
227 | ||||
228 | if (!isMIP && lpsolver == 2) | |||
229 | { | |||
230 | glp_iptcp iptcp; | |||
231 | glp_init_iptcp (&iptcp); | |||
232 | iptcp.msg_lev = par->msglev; | |||
233 | errnum = glp_interior (lp, &iptcp); | |||
234 | } | |||
235 | ||||
236 | if (errnum == 0) | |||
237 | { | |||
238 | if (isMIP) | |||
239 | { | |||
240 | *status = glp_mip_status (lp); | |||
241 | *fmin = glp_mip_obj_val (lp); | |||
242 | } | |||
243 | else | |||
244 | { | |||
245 | if (lpsolver == 1) | |||
246 | { | |||
247 | *status = glp_get_status (lp); | |||
248 | *fmin = glp_get_obj_val (lp); | |||
249 | } | |||
250 | else | |||
251 | { | |||
252 | *status = glp_ipt_status (lp); | |||
253 | *fmin = glp_ipt_obj_val (lp); | |||
254 | } | |||
255 | } | |||
256 | ||||
257 | if (isMIP) | |||
258 | { | |||
259 | for (int i = 0; i < n; i++) | |||
260 | xmin[i] = glp_mip_col_val (lp, i+1); | |||
261 | } | |||
262 | else | |||
263 | { | |||
264 | /* Primal values */ | |||
265 | for (int i = 0; i < n; i++) | |||
266 | { | |||
267 | if (lpsolver == 1) | |||
268 | xmin[i] = glp_get_col_prim (lp, i+1); | |||
269 | else | |||
270 | xmin[i] = glp_ipt_col_prim (lp, i+1); | |||
271 | } | |||
272 | ||||
273 | /* Dual values */ | |||
274 | for (int i = 0; i < m; i++) | |||
275 | { | |||
276 | if (lpsolver == 1) | |||
277 | lambda[i] = glp_get_row_dual (lp, i+1); | |||
278 | else | |||
279 | lambda[i] = glp_ipt_row_dual (lp, i+1); | |||
280 | } | |||
281 | ||||
282 | /* Reduced costs */ | |||
283 | for (int i = 0; i < glp_get_num_cols (lp); i++) | |||
284 | { | |||
285 | if (lpsolver == 1) | |||
286 | redcosts[i] = glp_get_col_dual (lp, i+1); | |||
287 | else | |||
288 | redcosts[i] = glp_ipt_col_dual (lp, i+1); | |||
289 | } | |||
290 | } | |||
291 | ||||
292 | *time = (clock () - t_start) / CLOCKS_PER_SEC1000000l; | |||
293 | } | |||
294 | ||||
295 | glp_delete_prob (lp); | |||
296 | ||||
297 | return errnum; | |||
298 | } | |||
299 | ||||
300 | #endif | |||
301 | ||||
302 | #define OCTAVE_GLPK_GET_REAL_PARAM(NAME, VAL)do { octave_value tmp = PARAM.getfield (NAME); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { VAL = tmp.scalar_value (); if (error_state) { error ("glpk: invalid value in PARAM." NAME) ; return retval; } } else { error ("glpk: invalid value in PARAM." NAME); return retval; } } } while (0) \ | |||
303 | do \ | |||
304 | { \ | |||
305 | octave_value tmp = PARAM.getfield (NAME); \ | |||
306 | \ | |||
307 | if (tmp.is_defined ()) \ | |||
308 | { \ | |||
309 | if (! tmp.is_empty ()) \ | |||
310 | { \ | |||
311 | VAL = tmp.scalar_value (); \ | |||
312 | \ | |||
313 | if (error_state) \ | |||
314 | { \ | |||
315 | error ("glpk: invalid value in PARAM." NAME); \ | |||
316 | return retval; \ | |||
317 | } \ | |||
318 | } \ | |||
319 | else \ | |||
320 | { \ | |||
321 | error ("glpk: invalid value in PARAM." NAME); \ | |||
322 | return retval; \ | |||
323 | } \ | |||
324 | } \ | |||
325 | } \ | |||
326 | while (0) | |||
327 | ||||
328 | #define OCTAVE_GLPK_GET_INT_PARAM(NAME, VAL)do { octave_value tmp = PARAM.getfield (NAME); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { VAL = tmp.int_value (); if (error_state ) { error ("glpk: invalid value in PARAM." NAME); return retval ; } } else { error ("glpk: invalid value in PARAM." NAME); return retval; } } } while (0) \ | |||
329 | do \ | |||
330 | { \ | |||
331 | octave_value tmp = PARAM.getfield (NAME); \ | |||
332 | \ | |||
333 | if (tmp.is_defined ()) \ | |||
334 | { \ | |||
335 | if (! tmp.is_empty ()) \ | |||
336 | { \ | |||
337 | VAL = tmp.int_value (); \ | |||
338 | \ | |||
339 | if (error_state) \ | |||
340 | { \ | |||
341 | error ("glpk: invalid value in PARAM." NAME); \ | |||
342 | return retval; \ | |||
343 | } \ | |||
344 | } \ | |||
345 | else \ | |||
346 | { \ | |||
347 | error ("glpk: invalid value in PARAM." NAME); \ | |||
348 | return retval; \ | |||
349 | } \ | |||
350 | } \ | |||
351 | } \ | |||
352 | while (0) | |||
353 | ||||
354 | DEFUN_DLD (__glpk__, args, ,octave_value_list F__glpk__ (const octave_value_list& args , int ); extern "C" octave_function * G__glpk__ (const octave_shlib & shl, bool relative) { octave_function *retval = 0; check_version ("api-v48+", "__glpk__"); if (! error_state) { octave_dld_function *fcn = octave_dld_function::create (F__glpk__, shl, "__glpk__" , "-*- texinfo -*-\n@deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\nUndocumented internal function.\n@end deftypefn" ); if (relative) fcn->mark_relative (); retval = fcn; } return retval; } octave_value_list F__glpk__ (const octave_value_list & args, int ) | |||
355 | "-*- texinfo -*-\n\octave_value_list F__glpk__ (const octave_value_list& args , int ); extern "C" octave_function * G__glpk__ (const octave_shlib & shl, bool relative) { octave_function *retval = 0; check_version ("api-v48+", "__glpk__"); if (! error_state) { octave_dld_function *fcn = octave_dld_function::create (F__glpk__, shl, "__glpk__" , "-*- texinfo -*-\n@deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\nUndocumented internal function.\n@end deftypefn" ); if (relative) fcn->mark_relative (); retval = fcn; } return retval; } octave_value_list F__glpk__ (const octave_value_list & args, int ) | |||
356 | @deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\n\octave_value_list F__glpk__ (const octave_value_list& args , int ); extern "C" octave_function * G__glpk__ (const octave_shlib & shl, bool relative) { octave_function *retval = 0; check_version ("api-v48+", "__glpk__"); if (! error_state) { octave_dld_function *fcn = octave_dld_function::create (F__glpk__, shl, "__glpk__" , "-*- texinfo -*-\n@deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\nUndocumented internal function.\n@end deftypefn" ); if (relative) fcn->mark_relative (); retval = fcn; } return retval; } octave_value_list F__glpk__ (const octave_value_list & args, int ) | |||
357 | Undocumented internal function.\n\octave_value_list F__glpk__ (const octave_value_list& args , int ); extern "C" octave_function * G__glpk__ (const octave_shlib & shl, bool relative) { octave_function *retval = 0; check_version ("api-v48+", "__glpk__"); if (! error_state) { octave_dld_function *fcn = octave_dld_function::create (F__glpk__, shl, "__glpk__" , "-*- texinfo -*-\n@deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\nUndocumented internal function.\n@end deftypefn" ); if (relative) fcn->mark_relative (); retval = fcn; } return retval; } octave_value_list F__glpk__ (const octave_value_list & args, int ) | |||
358 | @end deftypefn")octave_value_list F__glpk__ (const octave_value_list& args , int ); extern "C" octave_function * G__glpk__ (const octave_shlib & shl, bool relative) { octave_function *retval = 0; check_version ("api-v48+", "__glpk__"); if (! error_state) { octave_dld_function *fcn = octave_dld_function::create (F__glpk__, shl, "__glpk__" , "-*- texinfo -*-\n@deftypefn {Loadable Function} {[@var{values}] =} __glpk__ (@var{args})\nUndocumented internal function.\n@end deftypefn" ); if (relative) fcn->mark_relative (); retval = fcn; } return retval; } octave_value_list F__glpk__ (const octave_value_list & args, int ) | |||
359 | { | |||
360 | // The list of values to return. See the declaration in oct-obj.h | |||
361 | octave_value_list retval; | |||
362 | ||||
363 | #if defined (HAVE_GLPK1) | |||
364 | ||||
365 | int nrhs = args.length (); | |||
366 | ||||
367 | if (nrhs != 9) | |||
| ||||
368 | { | |||
369 | print_usage (); | |||
370 | return retval; | |||
371 | } | |||
372 | ||||
373 | //-- 1nd Input. A column array containing the objective function | |||
374 | //-- coefficients. | |||
375 | volatile int mrowsc = args(0).rows (); | |||
376 | ||||
377 | Matrix C (args(0).matrix_value ()); | |||
378 | ||||
379 | if (error_state) | |||
380 | { | |||
381 | error ("__glpk__: invalid value of C"); | |||
382 | return retval; | |||
383 | } | |||
384 | ||||
385 | double *c = C.fortran_vec (); | |||
386 | Array<int> rn; | |||
387 | Array<int> cn; | |||
388 | ColumnVector a; | |||
389 | volatile int mrowsA; | |||
390 | volatile int nz = 0; | |||
391 | ||||
392 | //-- 2nd Input. A matrix containing the constraints coefficients. | |||
393 | // If matrix A is NOT a sparse matrix | |||
394 | if (args(1).is_sparse_type ()) | |||
395 | { | |||
396 | SparseMatrix A = args(1).sparse_matrix_value (); // get the sparse matrix | |||
397 | ||||
398 | if (error_state) | |||
399 | { | |||
400 | error ("__glpk__: invalid value of A"); | |||
401 | return retval; | |||
402 | } | |||
403 | ||||
404 | mrowsA = A.rows (); | |||
405 | octave_idx_type Anc = A.cols (); | |||
406 | octave_idx_type Anz = A.nnz (); | |||
407 | rn.resize (dim_vector (Anz+1, 1)); | |||
408 | cn.resize (dim_vector (Anz+1, 1)); | |||
409 | a.resize (Anz+1, 0.0); | |||
410 | ||||
411 | if (Anc != mrowsc) | |||
412 | { | |||
413 | error ("__glpk__: invalid value of A"); | |||
414 | return retval; | |||
415 | } | |||
416 | ||||
417 | for (octave_idx_type j = 0; j < Anc; j++) | |||
418 | for (octave_idx_type i = A.cidx (j); i < A.cidx (j+1); i++) | |||
419 | { | |||
420 | nz++; | |||
421 | rn(nz) = A.ridx (i) + 1; | |||
422 | cn(nz) = j + 1; | |||
423 | a(nz) = A.data(i); | |||
424 | } | |||
425 | } | |||
426 | else | |||
427 | { | |||
428 | Matrix A (args(1).matrix_value ()); // get the matrix | |||
429 | ||||
430 | if (error_state) | |||
431 | { | |||
432 | error ("__glpk__: invalid value of A"); | |||
433 | return retval; | |||
434 | } | |||
435 | ||||
436 | mrowsA = A.rows (); | |||
437 | rn.resize (dim_vector (mrowsA*mrowsc+1, 1)); | |||
438 | cn.resize (dim_vector (mrowsA*mrowsc+1, 1)); | |||
439 | a.resize (mrowsA*mrowsc+1, 0.0); | |||
440 | ||||
441 | for (int i = 0; i < mrowsA; i++) | |||
442 | { | |||
443 | for (int j = 0; j < mrowsc; j++) | |||
444 | { | |||
445 | if (A(i,j) != 0) | |||
446 | { | |||
447 | nz++; | |||
448 | rn(nz) = i + 1; | |||
449 | cn(nz) = j + 1; | |||
450 | a(nz) = A(i,j); | |||
451 | } | |||
452 | } | |||
453 | } | |||
454 | ||||
455 | } | |||
456 | ||||
457 | //-- 3rd Input. A column array containing the right-hand side value | |||
458 | // for each constraint in the constraint matrix. | |||
459 | Matrix B (args(2).matrix_value ()); | |||
460 | ||||
461 | if (error_state) | |||
462 | { | |||
463 | error ("__glpk__: invalid value of B"); | |||
464 | return retval; | |||
465 | } | |||
466 | ||||
467 | double *b = B.fortran_vec (); | |||
468 | ||||
469 | //-- 4th Input. An array of length mrowsc containing the lower | |||
470 | //-- bound on each of the variables. | |||
471 | Matrix LB (args(3).matrix_value ()); | |||
472 | ||||
473 | if (error_state || LB.length () < mrowsc) | |||
474 | { | |||
475 | error ("__glpk__: invalid value of LB"); | |||
476 | return retval; | |||
477 | } | |||
478 | ||||
479 | double *lb = LB.fortran_vec (); | |||
480 | ||||
481 | //-- LB argument, default: Free | |||
482 | Array<int> freeLB (dim_vector (mrowsc, 1)); | |||
483 | for (int i = 0; i < mrowsc; i++) | |||
484 | { | |||
485 | if (xisinf (lb[i])) | |||
486 | { | |||
487 | freeLB(i) = 1; | |||
488 | lb[i] = -octave_Inf; | |||
489 | } | |||
490 | else | |||
491 | freeLB(i) = 0; | |||
492 | } | |||
493 | ||||
494 | //-- 5th Input. An array of at least length numcols containing the upper | |||
495 | //-- bound on each of the variables. | |||
496 | Matrix UB (args(4).matrix_value ()); | |||
497 | ||||
498 | if (error_state || UB.length () < mrowsc) | |||
499 | { | |||
500 | error ("__glpk__: invalid value of UB"); | |||
501 | return retval; | |||
502 | } | |||
503 | ||||
504 | double *ub = UB.fortran_vec (); | |||
505 | ||||
506 | Array<int> freeUB (dim_vector (mrowsc, 1)); | |||
507 | for (int i = 0; i < mrowsc; i++) | |||
508 | { | |||
509 | if (xisinf (ub[i])) | |||
510 | { | |||
511 | freeUB(i) = 1; | |||
512 | ub[i] = octave_Inf; | |||
513 | } | |||
514 | else | |||
515 | freeUB(i) = 0; | |||
516 | } | |||
517 | ||||
518 | //-- 6th Input. A column array containing the sense of each constraint | |||
519 | //-- in the constraint matrix. | |||
520 | charMatrix CTYPE (args(5).char_matrix_value ()); | |||
521 | ||||
522 | if (error_state) | |||
523 | { | |||
524 | error ("__glpk__: invalid value of CTYPE"); | |||
525 | return retval; | |||
526 | } | |||
527 | ||||
528 | char *ctype = CTYPE.fortran_vec (); | |||
529 | ||||
530 | //-- 7th Input. A column array containing the types of the variables. | |||
531 | charMatrix VTYPE (args(6).char_matrix_value ()); | |||
532 | ||||
533 | if (error_state) | |||
534 | { | |||
535 | error ("__glpk__: invalid value of VARTYPE"); | |||
536 | return retval; | |||
537 | } | |||
538 | ||||
539 | Array<int> vartype (dim_vector (mrowsc, 1)); | |||
540 | volatile int isMIP = 0; | |||
541 | for (int i = 0; i < mrowsc ; i++) | |||
542 | { | |||
543 | if (VTYPE(i,0) == 'I') | |||
544 | { | |||
545 | isMIP = 1; | |||
546 | vartype(i) = GLP_IV2; | |||
547 | } | |||
548 | else | |||
549 | vartype(i) = GLP_CV1; | |||
550 | } | |||
551 | ||||
552 | //-- 8th Input. Sense of optimization. | |||
553 | volatile int sense; | |||
554 | double SENSE = args(7).scalar_value (); | |||
555 | ||||
556 | if (error_state) | |||
557 | { | |||
558 | error ("__glpk__: invalid value of SENSE"); | |||
559 | return retval; | |||
560 | } | |||
561 | ||||
562 | if (SENSE >= 0) | |||
563 | sense = 1; | |||
564 | else | |||
565 | sense = -1; | |||
566 | ||||
567 | //-- 9th Input. A structure containing the control parameters. | |||
568 | octave_scalar_map PARAM = args(8).scalar_map_value (); | |||
569 | ||||
570 | if (error_state) | |||
571 | { | |||
572 | error ("__glpk__: invalid value of PARAM"); | |||
573 | return retval; | |||
574 | } | |||
575 | ||||
576 | control_params par; | |||
577 | ||||
578 | //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |||
579 | //-- Integer parameters | |||
580 | //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |||
581 | ||||
582 | //-- Level of messages output by the solver | |||
583 | par.msglev = 1; | |||
584 | OCTAVE_GLPK_GET_INT_PARAM ("msglev", par.msglev)do { octave_value tmp = PARAM.getfield ("msglev"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.msglev = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "msglev" ); return retval; } } else { error ("glpk: invalid value in PARAM." "msglev"); return retval; } } } while (0); | |||
585 | if (par.msglev < 0 || par.msglev > 3) | |||
586 | { | |||
587 | error ("__glpk__: PARAM.msglev must be 0 (no output) or 1 (error and warning messages only [default]) or 2 (normal output) or 3 (full output)"); | |||
588 | return retval; | |||
589 | } | |||
590 | ||||
591 | //-- scaling option | |||
592 | volatile int scale = 16; | |||
593 | OCTAVE_GLPK_GET_INT_PARAM ("scale", scale)do { octave_value tmp = PARAM.getfield ("scale"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { scale = tmp.int_value (); if ( error_state) { error ("glpk: invalid value in PARAM." "scale" ); return retval; } } else { error ("glpk: invalid value in PARAM." "scale"); return retval; } } } while (0); | |||
594 | if (scale < 0 || scale > 128) | |||
595 | { | |||
596 | error ("__glpk__: PARAM.scale must either be 128 (automatic selection of scaling options), or a bitwise or of: 1 (geometric mean scaling), 16 (equilibration scaling), 32 (round scale factors to power of two), 64 (skip if problem is well scaled"); | |||
597 | return retval; | |||
598 | } | |||
599 | ||||
600 | //-- Dual simplex option | |||
601 | par.dual = 1; | |||
602 | OCTAVE_GLPK_GET_INT_PARAM ("dual", par.dual)do { octave_value tmp = PARAM.getfield ("dual"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.dual = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "dual" ); return retval; } } else { error ("glpk: invalid value in PARAM." "dual"); return retval; } } } while (0); | |||
603 | if (par.dual < 1 || par.dual > 3) | |||
604 | { | |||
605 | error ("__glpk__: PARAM.dual must be 1 (use two-phase primal simplex [default]) or 2 (use two-phase dual simplex) or 3 (use two-phase dual simplex, and if it fails, switch to the primal simplex)"); | |||
606 | return retval; | |||
607 | } | |||
608 | ||||
609 | //-- Pricing option | |||
610 | par.price = 34; | |||
611 | OCTAVE_GLPK_GET_INT_PARAM ("price", par.price)do { octave_value tmp = PARAM.getfield ("price"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.price = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "price" ); return retval; } } else { error ("glpk: invalid value in PARAM." "price"); return retval; } } } while (0); | |||
612 | if (par.price != 17 && par.price != 34) | |||
613 | { | |||
614 | error ("__glpk__: PARAM.price must be 17 (textbook pricing) or 34 (steepest edge pricing [default])"); | |||
615 | return retval; | |||
616 | } | |||
617 | ||||
618 | //-- Simplex iterations limit | |||
619 | par.itlim = std::numeric_limits<int>::max (); | |||
620 | OCTAVE_GLPK_GET_INT_PARAM ("itlim", par.itlim)do { octave_value tmp = PARAM.getfield ("itlim"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.itlim = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "itlim" ); return retval; } } else { error ("glpk: invalid value in PARAM." "itlim"); return retval; } } } while (0); | |||
621 | ||||
622 | //-- Output frequency, in iterations | |||
623 | par.outfrq = 200; | |||
624 | OCTAVE_GLPK_GET_INT_PARAM ("outfrq", par.outfrq)do { octave_value tmp = PARAM.getfield ("outfrq"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.outfrq = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "outfrq" ); return retval; } } else { error ("glpk: invalid value in PARAM." "outfrq"); return retval; } } } while (0); | |||
625 | ||||
626 | //-- Branching heuristic option | |||
627 | par.branch = 4; | |||
628 | OCTAVE_GLPK_GET_INT_PARAM ("branch", par.branch)do { octave_value tmp = PARAM.getfield ("branch"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.branch = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "branch" ); return retval; } } else { error ("glpk: invalid value in PARAM." "branch"); return retval; } } } while (0); | |||
629 | if (par.branch < 1 || par.branch > 5) | |||
630 | { | |||
631 | error ("__glpk__: PARAM.branch must be 1 (first fractional variable) or 2 (last fractional variable) or 3 (most fractional variable) or 4 (heuristic by Driebeck and Tomlin [default]) or 5 (hybrid pseudocost heuristic)"); | |||
632 | return retval; | |||
633 | } | |||
634 | ||||
635 | //-- Backtracking heuristic option | |||
636 | par.btrack = 4; | |||
637 | OCTAVE_GLPK_GET_INT_PARAM ("btrack", par.btrack)do { octave_value tmp = PARAM.getfield ("btrack"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.btrack = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "btrack" ); return retval; } } else { error ("glpk: invalid value in PARAM." "btrack"); return retval; } } } while (0); | |||
638 | if (par.btrack < 1 || par.btrack > 4) | |||
639 | { | |||
640 | error ("__glpk__: PARAM.btrack must be 1 (depth first search) or 2 (breadth first search) or 3 (best local bound) or 4 (best projection heuristic [default]"); | |||
641 | return retval; | |||
642 | } | |||
643 | ||||
644 | //-- Presolver option | |||
645 | par.presol = 1; | |||
646 | OCTAVE_GLPK_GET_INT_PARAM ("presol", par.presol)do { octave_value tmp = PARAM.getfield ("presol"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.presol = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "presol" ); return retval; } } else { error ("glpk: invalid value in PARAM." "presol"); return retval; } } } while (0); | |||
647 | if (par.presol < 0 || par.presol > 1) | |||
648 | { | |||
649 | error ("__glpk__: PARAM.presol must be 0 (do NOT use LP presolver) or 1 (use LP presolver [default])"); | |||
650 | return retval; | |||
651 | } | |||
652 | ||||
653 | //-- LPsolver option | |||
654 | volatile int lpsolver = 1; | |||
655 | OCTAVE_GLPK_GET_INT_PARAM ("lpsolver", lpsolver)do { octave_value tmp = PARAM.getfield ("lpsolver"); if (tmp. is_defined ()) { if (! tmp.is_empty ()) { lpsolver = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "lpsolver"); return retval; } } else { error ("glpk: invalid value in PARAM." "lpsolver"); return retval; } } } while (0); | |||
656 | if (lpsolver < 1 || lpsolver > 2) | |||
657 | { | |||
658 | error ("__glpk__: PARAM.lpsolver must be 1 (simplex method) or 2 (interior point method)"); | |||
659 | return retval; | |||
660 | } | |||
661 | ||||
662 | //-- Ratio test option | |||
663 | par.rtest = 34; | |||
664 | OCTAVE_GLPK_GET_INT_PARAM ("rtest", par.rtest)do { octave_value tmp = PARAM.getfield ("rtest"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.rtest = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "rtest" ); return retval; } } else { error ("glpk: invalid value in PARAM." "rtest"); return retval; } } } while (0); | |||
665 | if (par.rtest != 17 && par.rtest != 34) | |||
666 | { | |||
667 | error ("__glpk__: PARAM.rtest must be 17 (standard ratio test) or 34 (Harris' two-pass ratio test [default])"); | |||
668 | return retval; | |||
669 | } | |||
670 | ||||
671 | par.tmlim = std::numeric_limits<int>::max (); | |||
672 | OCTAVE_GLPK_GET_INT_PARAM ("tmlim", par.tmlim)do { octave_value tmp = PARAM.getfield ("tmlim"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.tmlim = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "tmlim" ); return retval; } } else { error ("glpk: invalid value in PARAM." "tmlim"); return retval; } } } while (0); | |||
673 | ||||
674 | par.outdly = 0; | |||
675 | OCTAVE_GLPK_GET_INT_PARAM ("outdly", par.outdly)do { octave_value tmp = PARAM.getfield ("outdly"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.outdly = tmp.int_value () ; if (error_state) { error ("glpk: invalid value in PARAM." "outdly" ); return retval; } } else { error ("glpk: invalid value in PARAM." "outdly"); return retval; } } } while (0); | |||
676 | ||||
677 | //-- Save option | |||
678 | volatile int save_pb = 0; | |||
679 | OCTAVE_GLPK_GET_INT_PARAM ("save", save_pb)do { octave_value tmp = PARAM.getfield ("save"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { save_pb = tmp.int_value (); if (error_state) { error ("glpk: invalid value in PARAM." "save" ); return retval; } } else { error ("glpk: invalid value in PARAM." "save"); return retval; } } } while (0); | |||
680 | save_pb = save_pb != 0; | |||
681 | ||||
682 | //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |||
683 | //-- Real parameters | |||
684 | //-- ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ | |||
685 | ||||
686 | //-- Relative tolerance used to check if the current basic solution | |||
687 | //-- is primal feasible | |||
688 | par.tolbnd = 1e-7; | |||
689 | OCTAVE_GLPK_GET_REAL_PARAM ("tolbnd", par.tolbnd)do { octave_value tmp = PARAM.getfield ("tolbnd"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.tolbnd = tmp.scalar_value (); if (error_state) { error ("glpk: invalid value in PARAM." "tolbnd"); return retval; } } else { error ("glpk: invalid value in PARAM." "tolbnd"); return retval; } } } while (0); | |||
690 | ||||
691 | //-- Absolute tolerance used to check if the current basic solution | |||
692 | //-- is dual feasible | |||
693 | par.toldj = 1e-7; | |||
694 | OCTAVE_GLPK_GET_REAL_PARAM ("toldj", par.toldj)do { octave_value tmp = PARAM.getfield ("toldj"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.toldj = tmp.scalar_value ( ); if (error_state) { error ("glpk: invalid value in PARAM." "toldj" ); return retval; } } else { error ("glpk: invalid value in PARAM." "toldj"); return retval; } } } while (0); | |||
695 | ||||
696 | //-- Relative tolerance used to choose eligible pivotal elements of | |||
697 | //-- the simplex table in the ratio test | |||
698 | par.tolpiv = 1e-10; | |||
699 | OCTAVE_GLPK_GET_REAL_PARAM ("tolpiv", par.tolpiv)do { octave_value tmp = PARAM.getfield ("tolpiv"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.tolpiv = tmp.scalar_value (); if (error_state) { error ("glpk: invalid value in PARAM." "tolpiv"); return retval; } } else { error ("glpk: invalid value in PARAM." "tolpiv"); return retval; } } } while (0); | |||
700 | ||||
701 | par.objll = -std::numeric_limits<double>::max (); | |||
702 | OCTAVE_GLPK_GET_REAL_PARAM ("objll", par.objll)do { octave_value tmp = PARAM.getfield ("objll"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.objll = tmp.scalar_value ( ); if (error_state) { error ("glpk: invalid value in PARAM." "objll" ); return retval; } } else { error ("glpk: invalid value in PARAM." "objll"); return retval; } } } while (0); | |||
703 | ||||
704 | par.objul = std::numeric_limits<double>::max (); | |||
705 | OCTAVE_GLPK_GET_REAL_PARAM ("objul", par.objul)do { octave_value tmp = PARAM.getfield ("objul"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.objul = tmp.scalar_value ( ); if (error_state) { error ("glpk: invalid value in PARAM." "objul" ); return retval; } } else { error ("glpk: invalid value in PARAM." "objul"); return retval; } } } while (0); | |||
706 | ||||
707 | par.tolint = 1e-5; | |||
708 | OCTAVE_GLPK_GET_REAL_PARAM ("tolint", par.tolint)do { octave_value tmp = PARAM.getfield ("tolint"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.tolint = tmp.scalar_value (); if (error_state) { error ("glpk: invalid value in PARAM." "tolint"); return retval; } } else { error ("glpk: invalid value in PARAM." "tolint"); return retval; } } } while (0); | |||
709 | ||||
710 | par.tolobj = 1e-7; | |||
711 | OCTAVE_GLPK_GET_REAL_PARAM ("tolobj", par.tolobj)do { octave_value tmp = PARAM.getfield ("tolobj"); if (tmp.is_defined ()) { if (! tmp.is_empty ()) { par.tolobj = tmp.scalar_value (); if (error_state) { error ("glpk: invalid value in PARAM." "tolobj"); return retval; } } else { error ("glpk: invalid value in PARAM." "tolobj"); return retval; } } } while (0); | |||
712 | ||||
713 | //-- Assign pointers to the output parameters | |||
714 | ColumnVector xmin (mrowsc, octave_NA); | |||
715 | double fmin = octave_NA; | |||
716 | ColumnVector lambda (mrowsA, octave_NA); | |||
717 | ColumnVector redcosts (mrowsc, octave_NA); | |||
718 | double time; | |||
719 | int status, errnum = 0; | |||
720 | ||||
721 | int jmpret = setjmp (mark)_setjmp (mark); | |||
722 | ||||
723 | if (jmpret == 0) | |||
724 | errnum = glpk (sense, mrowsc, mrowsA, c, nz, rn.fortran_vec (), | |||
725 | cn.fortran_vec (), a.fortran_vec (), b, ctype, | |||
726 | freeLB.fortran_vec (), lb, freeUB.fortran_vec (), ub, | |||
727 | vartype.fortran_vec (), isMIP, lpsolver, save_pb, scale, | |||
728 | &par, xmin.fortran_vec (), &fmin, &status, | |||
729 | lambda.fortran_vec (), redcosts.fortran_vec (), &time); | |||
730 | ||||
731 | octave_scalar_map extra; | |||
732 | ||||
733 | if (! isMIP) | |||
734 | { | |||
735 | extra.assign ("lambda", lambda); | |||
736 | extra.assign ("redcosts", redcosts); | |||
737 | } | |||
738 | ||||
739 | extra.assign ("time", time); | |||
| ||||
740 | extra.assign ("status", status); | |||
741 | ||||
742 | retval(3) = extra; | |||
743 | retval(2) = errnum; | |||
744 | retval(1) = fmin; | |||
745 | retval(0) = xmin; | |||
746 | ||||
747 | #else | |||
748 | ||||
749 | gripe_not_supported ("glpk"); | |||
750 | ||||
751 | #endif | |||
752 | ||||
753 | return retval; | |||
754 | } | |||
755 | ||||
756 | /* | |||
757 | ## No test needed for internal helper function. | |||
758 | %!assert (1) | |||
759 | */ |