1 |
/* ASCEND modelling environment |
2 |
Copyright (C) 2007 Carnegie Mellon University |
3 |
|
4 |
This program is free software; you can redistribute it and/or modify |
5 |
it under the terms of the GNU General Public License as published by |
6 |
the Free Software Foundation; either version 2, or (at your option) |
7 |
any later version. |
8 |
|
9 |
This program is distributed in the hope that it will be useful, |
10 |
but WITHOUT ANY WARRANTY; without even the implied warranty of |
11 |
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
12 |
GNU General Public License for more details. |
13 |
|
14 |
You should have received a copy of the GNU General Public License |
15 |
along with this program; if not, write to the Free Software |
16 |
Foundation, Inc., 59 Temple Place - Suite 330, |
17 |
Boston, MA 02111-1307, USA. |
18 |
*//** |
19 |
@file |
20 |
Connection of the IPOPT optimisation solver into ASCEND. |
21 |
|
22 |
THIS IS STILL VERY MUCH UNDER DEVELOPMENT AND INCOMPLETE. I'VE ACTUALLY |
23 |
ONLY JUST STARTED WRITING IT by starting with asc_tron.c and modifying. |
24 |
|
25 |
The IPOPT solver is documented at http://projects.coin-or.org/Ipopt/ |
26 |
*//* |
27 |
ASCEND wrapper for IPOPT originally by John Pye, Jun 2007. |
28 |
*/ |
29 |
|
30 |
#include <utilities/config.h> |
31 |
|
32 |
#ifndef ASC_WITH_IPOPT |
33 |
# error "ASC_WITH_IPOPT must be defined in order to build this." |
34 |
#endif |
35 |
|
36 |
#include <solver/solver.h> |
37 |
|
38 |
#include <system/calc.h> |
39 |
#include <system/relman.h> |
40 |
#include <system/slv_stdcalls.h> |
41 |
#include <system/block.h> |
42 |
|
43 |
#include <utilities/ascConfig.h> |
44 |
#include <utilities/ascPanic.h> |
45 |
#include <utilities/ascMalloc.h> |
46 |
#include <utilities/ascDynaLoad.h> |
47 |
#include <utilities/mem.h> |
48 |
#include <utilities/ascEnvVar.h> |
49 |
#include <general/tm_time.h> |
50 |
#include <general/env.h> |
51 |
|
52 |
#include <ipopt/IpStdCInterface.h> |
53 |
|
54 |
ASC_DLLSPEC SolverRegisterFn ipopt_register; |
55 |
|
56 |
/*------------------------------------------------------------------------------ |
57 |
DATA STRUCTURES AND FORWARD DEFS |
58 |
*/ |
59 |
|
60 |
/** |
61 |
Documentation of solver options for IPOPT is at |
62 |
http://www.coin-or.org/Ipopt/documentation/node1.html |
63 |
*/ |
64 |
enum{ |
65 |
IPOPT_PARAM_TOL |
66 |
,IPOPT_PARAM_MAX_ITER |
67 |
,IPOPT_PARAM_SAFEEVAL |
68 |
,IPOPT_PARAM_MU_STRATEGY |
69 |
,IPOPT_PARAMS |
70 |
}; |
71 |
|
72 |
#define SYS(s) ((IpoptSystem *)(s)) |
73 |
|
74 |
struct IpoptSystemStruct{ |
75 |
|
76 |
/* |
77 |
Problem definition |
78 |
*/ |
79 |
slv_system_t slv; /* slv_system_t back-link */ |
80 |
struct rel_relation *obj; /* Objective function: NULL = none */ |
81 |
struct rel_relation *old_obj;/* Objective function: NULL = none */ |
82 |
struct var_variable **vlist; /* Variable list (NULL terminated) */ |
83 |
struct rel_relation **rlist; /* Relation list (NULL terminated) */ |
84 |
var_filter_t vfilt; |
85 |
|
86 |
/* |
87 |
Solver information |
88 |
*/ |
89 |
int32 presolved; /* ? Has the system been presolved */ |
90 |
int32 resolve; /* ? Has the system been resolved */ |
91 |
slv_parameters_t p; /* Parameters */ |
92 |
slv_status_t s; /* Status (as of iteration end) */ |
93 |
|
94 |
int32 cap; /* Order of matrix/vectors */ |
95 |
int32 rank; /* Symbolic rank of problem */ |
96 |
int32 vused; /* Free and incident variables */ |
97 |
int32 vtot; /* length of varlist */ |
98 |
int32 rused; /* Included relations */ |
99 |
int32 rtot; /* length of rellist */ |
100 |
double clock; /* CPU time */ |
101 |
|
102 |
int32 calc_ok; |
103 |
double obj_val; |
104 |
|
105 |
void *parm_array[IPOPT_PARAMS]; |
106 |
struct slv_parameter pa[IPOPT_PARAMS]; |
107 |
|
108 |
/* |
109 |
IPOPT DATA |
110 |
*/ |
111 |
Index n; /* number of variables */ |
112 |
Index m; /* number of constraints */ |
113 |
|
114 |
int nnzJ; /* number of non zeros in the jacobian of the constraints */ |
115 |
int nnzH; /* number of non-zeros in the hessian of the objective */ |
116 |
|
117 |
Number* x_L; /* lower bounds on x */ |
118 |
Number* x_U; /* upper bounds on x */ |
119 |
Number* g_L; /* lower bounds on g */ |
120 |
Number* g_U; /* upper bounds on g */ |
121 |
|
122 |
IpoptProblem nlp; /* IpoptProblem */ |
123 |
|
124 |
enum ApplicationReturnStatus status; /* Solve return code */ |
125 |
Number* x; /* starting point and solution vector */ |
126 |
Number* mult_x_L; /* lower bound multipliers at the solution */ |
127 |
Number* mult_x_U; /* upper bound multipliers at the solution */ |
128 |
Index i; /* generic counter */ |
129 |
}; |
130 |
|
131 |
typedef struct IpoptSystemStruct IpoptSystem; |
132 |
|
133 |
static int ipopt_get_default_parameters(slv_system_t server, SlvClientToken asys, slv_parameters_t *parameters); |
134 |
|
135 |
static void ipopt_iteration_begins(IpoptSystem *sys); |
136 |
static void ipopt_iteration_ends(IpoptSystem *sys); |
137 |
|
138 |
/*------------------------------------------------------------------------------ |
139 |
SYSTEM SETUP/DESTROY AND SOLVER ELIGIBILITY |
140 |
*/ |
141 |
|
142 |
static SlvClientToken ipopt_create(slv_system_t server, int32*statusindex){ |
143 |
IpoptSystem *sys; |
144 |
|
145 |
sys = ASC_NEW_CLEAR(IpoptSystem); |
146 |
if(sys==NULL){ |
147 |
*statusindex = 1; |
148 |
return sys; |
149 |
} |
150 |
|
151 |
sys->p.parms = sys->pa; |
152 |
sys->p.dynamic_parms = 0; |
153 |
ipopt_get_default_parameters(server,(SlvClientToken)sys,&(sys->p)); |
154 |
sys->p.whose = (*statusindex); |
155 |
|
156 |
sys->presolved = 0; |
157 |
sys->resolve = 0; |
158 |
|
159 |
sys->n = -1; |
160 |
sys->m = -1; |
161 |
|
162 |
sys->s.ok = TRUE; |
163 |
sys->s.calc_ok = TRUE; |
164 |
sys->s.costsize = 0; |
165 |
sys->s.cost = NULL; /*redundant, but sanity preserving */ |
166 |
sys->s.block.number_of = 1; |
167 |
sys->s.block.current_block = 0; |
168 |
sys->s.block.current_reordered_block = 0; |
169 |
sys->s.block.current_size = 0; |
170 |
sys->s.block.previous_total_size = 0; |
171 |
sys->s.block.iteration = 0; |
172 |
sys->s.block.funcs = 0; |
173 |
sys->s.block.jacs = 0; |
174 |
sys->s.block.cpu_elapsed = 0; |
175 |
sys->s.block.functime = 0; |
176 |
sys->s.block.jactime = 0; |
177 |
sys->s.block.residual = 0; |
178 |
|
179 |
/** @TODO set up sys->vfilt */ |
180 |
|
181 |
sys->vlist = slv_get_solvers_var_list(server); |
182 |
sys->rlist = slv_get_solvers_rel_list(server); |
183 |
sys->obj = slv_get_obj_relation(server); |
184 |
if(sys->vlist == NULL) { |
185 |
ASC_FREE(sys); |
186 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"IPOPT called with no variables."); |
187 |
*statusindex = -2; |
188 |
return NULL; |
189 |
} |
190 |
if(sys->rlist == NULL && sys->obj == NULL) { |
191 |
ASC_FREE(sys); |
192 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"IPOPT called with no relations or objective."); |
193 |
*statusindex = -1; |
194 |
return NULL; |
195 |
} |
196 |
|
197 |
/* do nothing with the objective list, pars, bounds, extrels, etc etc */ |
198 |
|
199 |
slv_check_var_initialization(server); |
200 |
*statusindex = 0; |
201 |
return((SlvClientToken)sys); |
202 |
} |
203 |
|
204 |
static int32 ipopt_destroy(slv_system_t server, SlvClientToken asys){ |
205 |
UNUSED_PARAMETER(server); |
206 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"ipopt_destroy not implemented"); |
207 |
return 1; |
208 |
} |
209 |
|
210 |
static int32 ipopt_eligible_solver(slv_system_t server){ |
211 |
UNUSED_PARAMETER(server); |
212 |
|
213 |
/// TODO check that there is a MAXIMIZE or MINIMIZE statement |
214 |
/// TODO check that there are no discrete-valued variables |
215 |
/// TODO check that there are no WHENs or CONDITIONALs |
216 |
/// TODO check anything else? |
217 |
|
218 |
ERROR_REPORTER_HERE(ASC_PROG_WARNING,"ipopt_eligible_solver not implemented"); |
219 |
return 1; |
220 |
} |
221 |
|
222 |
/*------------------------------------------------------------------------------ |
223 |
SOLVER PARAMETERS |
224 |
*/ |
225 |
|
226 |
static |
227 |
int32 ipopt_get_default_parameters(slv_system_t server, SlvClientToken asys |
228 |
,slv_parameters_t *parameters |
229 |
){ |
230 |
IpoptSystem *sys = NULL; |
231 |
struct slv_parameter *new_parms = NULL; |
232 |
int32 make_macros = 0; |
233 |
|
234 |
if(server != NULL && asys != NULL) { |
235 |
sys = SYS(asys); |
236 |
make_macros = 1; |
237 |
} |
238 |
|
239 |
if(parameters->parms == NULL) { |
240 |
new_parms = ASC_NEW_ARRAY_OR_NULL(struct slv_parameter,IPOPT_PARAMS); |
241 |
if(new_parms == NULL) { |
242 |
return -1; |
243 |
} |
244 |
parameters->parms = new_parms; |
245 |
parameters->dynamic_parms = 1; |
246 |
} |
247 |
|
248 |
parameters->num_parms = 0; |
249 |
|
250 |
slv_param_int(parameters,IPOPT_PARAM_MAX_ITER |
251 |
,(SlvParameterInitInt){{"max_iter" |
252 |
,"Maximum number of iterations",2 |
253 |
,"The algorithm terminates with an error message if the number of iterations exceeded this number." |
254 |
}, 3000, 0, 100000000} |
255 |
); |
256 |
|
257 |
slv_param_real(parameters,IPOPT_PARAM_TOL |
258 |
,(SlvParameterInitReal){{"tol" |
259 |
,"Desired convergence tolerance (relative)",2 |
260 |
,"Determines the convergence tolerance for the algorithm. The algorithm" |
261 |
" terminates successfully, if the (scaled) NLP error becomes smaller" |
262 |
" than this value, and if the (absolute) criteria according to " |
263 |
" 'dual_inf_tol', 'primal_inf_tol', and 'cmpl_inf_tol' are met. (This" |
264 |
" is epsilon_tol in Eqn. (6) in implementation paper). See also " |
265 |
" 'acceptable_tol' as a second termination criterion. Note, some other" |
266 |
" algorithmic features also use this quantity to determine thresholds" |
267 |
" etc." |
268 |
}, 1e-8, 0, 1e20} |
269 |
); |
270 |
|
271 |
slv_param_char(parameters,IPOPT_PARAM_MU_STRATEGY |
272 |
,(SlvParameterInitChar){{"mu_strategy" |
273 |
,"Update strategy for barrier parameter",6 |
274 |
,"Determines which barrier parameter update strategy is to be used." |
275 |
" 'monotone' is the monotone (Fiacco-McCormick) strategy;" |
276 |
" 'adaptive' is the adaptive update strategy." |
277 |
}, "monotone"}, (char *[]){ |
278 |
"monotone","adaptive",NULL |
279 |
} |
280 |
); |
281 |
|
282 |
slv_param_bool(parameters,IPOPT_PARAM_SAFEEVAL |
283 |
,(SlvParameterInitBool){{"safeeval" |
284 |
,"Use safe evaluation?",1 |
285 |
,"Use 'safe' function evaluation routines (TRUE) or allow ASCEND to " |
286 |
"throw SIGFPE errors which will then halt integration (FALSE)." |
287 |
}, FALSE} |
288 |
); |
289 |
|
290 |
asc_assert(parameters->num_parms==IPOPT_PARAMS); |
291 |
|
292 |
return 1; |
293 |
} |
294 |
|
295 |
static void ipopt_get_parameters(slv_system_t server, SlvClientToken asys |
296 |
, slv_parameters_t *parameters |
297 |
){ |
298 |
IpoptSystem *sys; |
299 |
UNUSED_PARAMETER(server); |
300 |
|
301 |
sys = SYS(asys); |
302 |
//if(check_system(sys)) return; |
303 |
mem_copy_cast(&(sys->p),parameters,sizeof(slv_parameters_t)); |
304 |
} |
305 |
|
306 |
|
307 |
static void ipopt_set_parameters(slv_system_t server, SlvClientToken asys |
308 |
,slv_parameters_t *parameters |
309 |
){ |
310 |
IpoptSystem *sys; |
311 |
UNUSED_PARAMETER(server); |
312 |
|
313 |
sys = SYS(asys); |
314 |
//if (check_system(sys)) return; |
315 |
mem_copy_cast(parameters,&(sys->p),sizeof(slv_parameters_t)); |
316 |
} |
317 |
|
318 |
/*------------------------------------------------------------------------------ |
319 |
EVALUATION AND RESULT HOOK FUNCTIONS |
320 |
*/ |
321 |
|
322 |
/** |
323 |
update the model with new 'x' vector. |
324 |
@return 0 on success. |
325 |
*/ |
326 |
int ipopt_update_model(IpoptSystem *sys, const double *x){ |
327 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Not implemented"); |
328 |
return 1; /* error */ |
329 |
} |
330 |
|
331 |
/** Function to evaluate the objective function f(x). |
332 |
@return 1 on success, 0 on failure |
333 |
|
334 |
@param n (in), the number of variables in the problem (dimension of 'x'). |
335 |
@param x (in), the values for the primal variables, 'x' , at which 'f(x)' is to be evaluated. |
336 |
@param new_x (in), false if any evaluation method was previously called with the same values in 'x', true otherwise. |
337 |
@param obj_value (out) the value of the objective function ('f(x)'). |
338 |
*/ |
339 |
Bool ipopt_eval_f(Index n, Number *x, Bool new_x, Number *obj_value, void *user_data){ |
340 |
IpoptSystem *sys; |
341 |
sys = SYS(user_data); |
342 |
int res; |
343 |
|
344 |
asc_assert(n==sys->n); |
345 |
asc_assert(sys->obj!=NULL); |
346 |
|
347 |
if(new_x){ |
348 |
res = ipopt_update_model(sys,x); |
349 |
if(res)return 0; /* fail model update */ |
350 |
} |
351 |
|
352 |
sys->calc_ok = TRUE; |
353 |
|
354 |
*obj_value = relman_eval(sys->obj,&(sys->calc_ok),SLV_PARAM_BOOL(&(sys->p),IPOPT_PARAM_SAFEEVAL)); |
355 |
|
356 |
return sys->calc_ok; |
357 |
} |
358 |
|
359 |
/** |
360 |
@return 1 on success |
361 |
*/ |
362 |
Bool ipopt_eval_grad_f(Index n, Number* x, Bool new_x, Number* grad_f, void *user_data){ |
363 |
IpoptSystem *sys; |
364 |
sys = SYS(user_data); |
365 |
int j, res, len; |
366 |
int count; |
367 |
double *derivatives; |
368 |
int *variables; |
369 |
static var_filter_t vfilter = { |
370 |
VAR_ACTIVE | VAR_INCIDENT | VAR_SVAR |
371 |
,VAR_ACTIVE | VAR_INCIDENT | VAR_SVAR | VAR_FIXED |
372 |
}; |
373 |
|
374 |
asc_assert(n==sys->n); |
375 |
asc_assert(sys->obj); |
376 |
|
377 |
if(new_x){ |
378 |
res = ipopt_update_model(sys,x); |
379 |
if(res)return 0; /* fail model update */ |
380 |
} |
381 |
|
382 |
|
383 |
/* evaluate grad_f(x) somehow */ |
384 |
for(j=0; j<n; ++j){ |
385 |
grad_f[j] = 0; |
386 |
} |
387 |
|
388 |
len = rel_n_incidences(sys->obj); |
389 |
variables = ASC_NEW_ARRAY(int,len); |
390 |
derivatives = ASC_NEW_ARRAY(double,len); |
391 |
|
392 |
relman_diff2( |
393 |
sys->obj,&vfilter,derivatives,variables |
394 |
, &count,SLV_PARAM_BOOL(&(sys->p),IPOPT_PARAM_SAFEEVAL) |
395 |
); |
396 |
|
397 |
for(j=0; j<len; ++j){ |
398 |
grad_f[variables[j]] = derivatives[j]; |
399 |
} |
400 |
|
401 |
ASC_FREE(variables); |
402 |
ASC_FREE(derivatives); |
403 |
|
404 |
return 1; /* success, presumably */ |
405 |
} |
406 |
|
407 |
Bool ipopt_eval_g(Index n, Number* x, Bool new_x, Index m, Number *g, void *user_data){ |
408 |
IpoptSystem *sys; |
409 |
sys = SYS(user_data); |
410 |
int i, res; |
411 |
|
412 |
asc_assert(n==sys->n); |
413 |
asc_assert(m==sys->m); |
414 |
|
415 |
if(new_x){ |
416 |
res = ipopt_update_model(sys,x); |
417 |
if(res)return 0; /* fail model update */ |
418 |
} |
419 |
|
420 |
for(i=0; i<m; ++i){ |
421 |
g[i] = 0; |
422 |
} |
423 |
|
424 |
return 0; /* fail: not yet implemented */ |
425 |
} |
426 |
|
427 |
Bool ipopt_eval_jac_g(Index n, Number* x, Bool new_x, Index m |
428 |
, Index nele_jac, Index* iRow, Index *jCol, Number* values |
429 |
, void *user_data |
430 |
){ |
431 |
IpoptSystem *sys; |
432 |
sys = SYS(user_data); |
433 |
int i,res; |
434 |
|
435 |
asc_assert(n==sys->n); |
436 |
asc_assert(nele_jac==sys->nnzJ); |
437 |
asc_assert(m==sys->m); |
438 |
|
439 |
if(new_x){ |
440 |
res = ipopt_update_model(sys,x); |
441 |
if(res)return 0; /* fail model update */ |
442 |
} |
443 |
|
444 |
/* loop through the constraints */ |
445 |
for(i=0; i<m; ++i){ |
446 |
/* get derivatives for constraint i */ |
447 |
/* insert the derivatives into the matrix in row i, columns j */ |
448 |
} |
449 |
|
450 |
return 0; /* fail: not yet implemented */ |
451 |
} |
452 |
|
453 |
Bool ipopt_eval_h(Index n, Number* x, Bool new_x |
454 |
, Number obj_factor, Index m, Number* lambda |
455 |
, Bool new_lambda, Index nele_hess, Index* iRow |
456 |
, Index* jCol, Number* values |
457 |
, void *user_data |
458 |
){ |
459 |
/* not sure about this one yet: do all the 2nd deriv things work for this? */ |
460 |
return 0; /* fail: not yet implemented */ |
461 |
} |
462 |
|
463 |
/*------------------------------------------------------------------------------ |
464 |
SOLVE ROUTINES |
465 |
*/ |
466 |
|
467 |
static int ipopt_solve(slv_system_t server, SlvClientToken asys){ |
468 |
IpoptSystem *sys; |
469 |
UNUSED_PARAMETER(server); |
470 |
enum ApplicationReturnStatus status; |
471 |
int ret, i, j; |
472 |
struct var_variable *var; |
473 |
|
474 |
sys = SYS(asys); |
475 |
|
476 |
double *x, *x_L, *x_U, *g_L, *g_U, *mult_x_L, *mult_x_U; |
477 |
|
478 |
/* set the number of variables and allocate space for the bounds */ |
479 |
x_L = ASC_NEW_ARRAY(Number,sys->n); |
480 |
x_U = ASC_NEW_ARRAY(Number,sys->n); |
481 |
|
482 |
/* @TODO set the values for the variable bounds */ |
483 |
int jj = 0; |
484 |
for(j = 0; j < sys->vtot; j++){ |
485 |
var = sys->vlist[j]; |
486 |
if(var_apply_filter(var,&(sys->vfilt))){ |
487 |
x_L[jj] = var_lower_bound(var); |
488 |
x_U[jj] = var_upper_bound(var); |
489 |
jj++; |
490 |
} |
491 |
} |
492 |
|
493 |
/** @TODO set bounds on the constraints? */ |
494 |
|
495 |
/* create the IpoptProblem */ |
496 |
sys->nlp = CreateIpoptProblem(sys->n, x_L, x_U, sys->m, g_L, g_U, sys->nnzJ, sys->nnzH, 0/*index style=C*/, |
497 |
&ipopt_eval_f, &ipopt_eval_g, &ipopt_eval_grad_f, |
498 |
&ipopt_eval_jac_g, &ipopt_eval_h |
499 |
); |
500 |
|
501 |
/* We can free the memory now - the values for the bounds have been |
502 |
copied internally in CreateIpoptProblem */ |
503 |
ASC_FREE(x_L); |
504 |
ASC_FREE(x_U); |
505 |
ASC_FREE(g_L); |
506 |
ASC_FREE(g_U); |
507 |
|
508 |
/* set some options */ |
509 |
AddIpoptNumOption(sys->nlp, "tol", SLV_PARAM_BOOL(&(sys->p),IPOPT_PARAM_TOL)); |
510 |
AddIpoptStrOption(sys->nlp, "mu_strategy", SLV_PARAM_CHAR(&(sys->p),IPOPT_PARAM_MU_STRATEGY)); |
511 |
|
512 |
/* initial values */ |
513 |
x = ASC_NEW_ARRAY(Number, sys->n); |
514 |
|
515 |
/** @TODO get values of 'x' from the model */ |
516 |
|
517 |
/* allocate space to store the bound multipliers at the solution */ |
518 |
mult_x_L = ASC_NEW_ARRAY(Number, sys->n); |
519 |
mult_x_U = ASC_NEW_ARRAY(Number, sys->n); |
520 |
|
521 |
/* solve the problem */ |
522 |
status = IpoptSolve(sys->nlp, x, NULL, &sys->obj_val, NULL, mult_x_L, mult_x_U, NULL); |
523 |
|
524 |
/** @TODO update the sys->s.xxxxx flags based on value of 'status' */ |
525 |
|
526 |
if (status == Solve_Succeeded) { |
527 |
CONSOLE_DEBUG("Solution of the primal variables, x"); |
528 |
for (i=0; i<sys->n; i++) { |
529 |
CONSOLE_DEBUG(" x[%d] = %e\n", i, x[i]); |
530 |
} |
531 |
|
532 |
CONSOLE_DEBUG("Solution of the bound multipliers, z_L and z_U"); |
533 |
for (i=0; i<sys->n; i++) { |
534 |
CONSOLE_DEBUG(" z_L[%d] = %e", i, mult_x_L[i]); |
535 |
} |
536 |
for (i=0; i<sys->n; i++) { |
537 |
CONSOLE_DEBUG(" z_U[%d] = %e", i, mult_x_U[i]); |
538 |
} |
539 |
|
540 |
CONSOLE_DEBUG("Objective value"); |
541 |
CONSOLE_DEBUG(" f(x*) = %e", sys->obj_val); |
542 |
|
543 |
ret = 0; /* success */ |
544 |
}else{ |
545 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Failed solve, unknown status"); |
546 |
ret = 1; /* failure */ |
547 |
} |
548 |
|
549 |
/* free allocated memory */ |
550 |
FreeIpoptProblem(sys->nlp); |
551 |
ASC_FREE(x); |
552 |
ASC_FREE(mult_x_L); |
553 |
ASC_FREE(mult_x_U); |
554 |
|
555 |
return ret; |
556 |
} |
557 |
|
558 |
static int ipopt_presolve(slv_system_t server, SlvClientToken asys){ |
559 |
IpoptSystem *sys; |
560 |
|
561 |
CONSOLE_DEBUG("PRESOLVE"); |
562 |
|
563 |
sys = SYS(asys); |
564 |
ipopt_iteration_begins(sys); |
565 |
//check_system(sys); |
566 |
|
567 |
asc_assert(sys->vlist && sys->rlist); |
568 |
|
569 |
sys->obj = slv_get_obj_relation(server); /*may have changed objective*/ |
570 |
if(!sys->obj){ |
571 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"No objective function was specified"); |
572 |
return -3; |
573 |
} |
574 |
|
575 |
#if 0 |
576 |
if(sys->presolved > 0) { /* system has been presolved before */ |
577 |
if(!conopt_dof_changed(sys) /*no changes in fixed or included flags*/ |
578 |
&& sys->p.partition == sys->J.old_partition |
579 |
&& sys->obj == sys->old_obj |
580 |
){ |
581 |
matrix_creation_needed = 0; |
582 |
CONOPT_CONSOLE_DEBUG("YOU JUST AVOIDED MATRIX DESTRUCTION/CREATION"); |
583 |
} |
584 |
} |
585 |
#endif |
586 |
|
587 |
#if 0 |
588 |
// check all this... |
589 |
|
590 |
/* set all relations as being 'unsatisfied' to start with... */ |
591 |
rp=sys->rlist; |
592 |
for( ind = 0; ind < sys->rtot; ++ind ) { |
593 |
rel_set_satisfied(rp[ind],FALSE); |
594 |
} |
595 |
|
596 |
if( matrix_creation_needed ) { |
597 |
|
598 |
cap = slv_get_num_solvers_rels(server); |
599 |
sys->cap = slv_get_num_solvers_vars(server); |
600 |
sys->cap = MAX(sys->cap,cap); |
601 |
vp=sys->vlist; |
602 |
for( ind = 0; ind < sys->vtot; ++ind ) { |
603 |
var_set_in_block(vp[ind],FALSE); |
604 |
} |
605 |
rp=sys->rlist; |
606 |
for( ind = 0; ind < sys->rtot; ++ind ) { |
607 |
rel_set_in_block(rp[ind],FALSE); |
608 |
/* rel_set_satisfied(rp[ind],FALSE); */ |
609 |
} |
610 |
|
611 |
sys->presolved = 1; /* full presolve recognized here */ |
612 |
sys->resolve = 0; /* initialize resolve flag */ |
613 |
|
614 |
/* @TODO calculate nnz for jacobian matrix of constraints */ |
615 |
|
616 |
/* provide sparsity structure for jacobian */ |
617 |
|
618 |
/** @TODO calculate nnz for hessian matrix */ |
619 |
|
620 |
/* provide sparsity structure for hessian matrix */ |
621 |
|
622 |
sys->J.old_partition = sys->p.partition; |
623 |
sys->old_obj = sys->obj; |
624 |
|
625 |
slv_sort_rels_and_vars(server,&(sys->con.m),&(sys->con.n)); |
626 |
CONOPT_CONSOLE_DEBUG("FOUND %d CONSTRAINTS AND %d VARS",sys->con.m,sys->con.n); |
627 |
if (sys->obj != NULL) { |
628 |
CONOPT_CONSOLE_DEBUG("ADDING OBJECT AS A ROW"); |
629 |
sys->con.m++; /* treat objective as a row */ |
630 |
} |
631 |
|
632 |
cntvect = ASC_NEW_ARRAY(int,COIDEF_Size()); |
633 |
COIDEF_Ini(cntvect); |
634 |
sys->con.cntvect = cntvect; |
635 |
CONOPT_CONSOLE_DEBUG("NUMBER OF CONSTRAINTS = %d",sys->con.m); |
636 |
COIDEF_NumVar(cntvect, &(sys->con.n)); |
637 |
COIDEF_NumCon(cntvect, &(sys->con.m)); |
638 |
sys->con.nz = num_jacobian_nonzeros(sys, &(sys->con.maxrow)); |
639 |
COIDEF_NumNZ(cntvect, &(sys->con.nz)); |
640 |
COIDEF_NumNlNz(cntvect, &(sys->con.nz)); |
641 |
|
642 |
sys->con.base = 0; |
643 |
COIDEF_Base(cntvect,&(sys->con.base)); |
644 |
COIDEF_ErrLim(cntvect, &(DOMLIM)); |
645 |
COIDEF_ItLim(cntvect, &(ITER_LIMIT)); |
646 |
|
647 |
if(sys->obj!=NULL){ |
648 |
sys->con.optdir = relman_obj_direction(sys->obj); |
649 |
sys->con.objcon = sys->con.m - 1; /* objective will be last row */ |
650 |
CONOPT_CONSOLE_DEBUG("SETTING OBJECTIVE CONSTRAINT TO BE %d",sys->con.objcon); |
651 |
}else{ |
652 |
sys->con.optdir = 0; |
653 |
sys->con.objcon = 0; |
654 |
} |
655 |
COIDEF_OptDir(cntvect, &(sys->con.optdir)); |
656 |
COIDEF_ObjCon(cntvect, &(sys->con.objcon)); |
657 |
|
658 |
temp = 0; |
659 |
COIDEF_StdOut(cntvect, &temp); |
660 |
|
661 |
COIDEF_ReadMatrix(cntvect, &conopt_readmatrix); |
662 |
COIDEF_FDEval(cntvect, &conopt_fdeval); |
663 |
COIDEF_Option(cntvect, &conopt_option); |
664 |
COIDEF_Solution(cntvect, &conopt_solution); |
665 |
COIDEF_Status(cntvect, &conopt_status); |
666 |
COIDEF_Message(cntvect, &asc_conopt_message); |
667 |
COIDEF_ErrMsg(cntvect, &conopt_errmsg); |
668 |
COIDEF_Progress(cntvect, &asc_conopt_progress); |
669 |
|
670 |
int debugfv = 1; |
671 |
COIDEF_DebugFV(cntvect, &debugfv); |
672 |
|
673 |
destroy_vectors(sys); |
674 |
destroy_matrices(sys); |
675 |
create_matrices(server,sys); |
676 |
create_vectors(sys); |
677 |
|
678 |
sys->s.block.current_reordered_block = -2; |
679 |
} |
680 |
|
681 |
//... |
682 |
|
683 |
if( matrix_creation_needed ) { |
684 |
destroy_array(sys->s.cost); |
685 |
sys->s.cost = create_zero_array(sys->s.costsize,struct slv_block_cost); |
686 |
for( ind = 0; ind < sys->s.costsize; ++ind ) { |
687 |
sys->s.cost[ind].reorder_method = -1; |
688 |
} |
689 |
} else { |
690 |
reset_cost(sys->s.cost,sys->s.costsize); |
691 |
} |
692 |
|
693 |
#endif |
694 |
|
695 |
/* Reset status */ |
696 |
sys->s.iteration = 0; |
697 |
sys->s.cpu_elapsed = 0.0; |
698 |
sys->s.converged = sys->s.diverged = sys->s.inconsistent = FALSE; |
699 |
sys->s.block.previous_total_size = 0; |
700 |
sys->s.costsize = 1+sys->s.block.number_of; |
701 |
|
702 |
/* set to go to first unconverged block */ |
703 |
sys->s.block.current_block = -1; |
704 |
sys->s.block.current_size = 0; |
705 |
sys->s.calc_ok = TRUE; |
706 |
sys->s.block.iteration = 0; |
707 |
sys->obj_val = MAXDOUBLE/2000.0; |
708 |
|
709 |
ipopt_iteration_ends(sys); |
710 |
sys->s.cost[sys->s.block.number_of].time=sys->s.cpu_elapsed; |
711 |
|
712 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"presolve completed"); |
713 |
return 0; |
714 |
} |
715 |
|
716 |
/** |
717 |
Prepare sys for entering an iteration, increasing the iteration counts |
718 |
and starting the clock. |
719 |
*/ |
720 |
static void ipopt_iteration_begins(IpoptSystem *sys){ |
721 |
sys->clock = tm_cpu_time(); |
722 |
++(sys->s.block.iteration); |
723 |
++(sys->s.iteration); |
724 |
} |
725 |
|
726 |
|
727 |
/* |
728 |
Prepare sys for exiting an iteration, stopping the clock and recording |
729 |
the cpu time. |
730 |
*/ |
731 |
static void ipopt_iteration_ends(IpoptSystem *sys){ |
732 |
double cpu_elapsed; /* elapsed this iteration */ |
733 |
|
734 |
cpu_elapsed = (double)(tm_cpu_time() - sys->clock); |
735 |
sys->s.block.cpu_elapsed += cpu_elapsed; |
736 |
sys->s.cpu_elapsed += cpu_elapsed; |
737 |
} |
738 |
|
739 |
|
740 |
|
741 |
static int ipopt_iterate(slv_system_t server, SlvClientToken asys){ |
742 |
UNUSED_PARAMETER(server); |
743 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Not implemented"); |
744 |
return 1; |
745 |
} |
746 |
|
747 |
static int ipopt_resolve(slv_system_t server, SlvClientToken asys){ |
748 |
UNUSED_PARAMETER(server); |
749 |
|
750 |
/* if implementing this, use the 'warm start' thing in IPOPT */ |
751 |
|
752 |
/* provide initial values of the 'multipliers' */ |
753 |
|
754 |
ERROR_REPORTER_HERE(ASC_PROG_ERR,"Not implemented"); |
755 |
return 1; |
756 |
} |
757 |
|
758 |
static const SlvFunctionsT ipopt_internals = { |
759 |
67 |
760 |
,"IPOPT" |
761 |
,ipopt_create |
762 |
,ipopt_destroy |
763 |
,ipopt_eligible_solver |
764 |
,ipopt_get_default_parameters |
765 |
,ipopt_get_parameters |
766 |
,ipopt_set_parameters |
767 |
,NULL |
768 |
,ipopt_solve |
769 |
,ipopt_presolve |
770 |
,ipopt_iterate |
771 |
,ipopt_resolve |
772 |
,NULL |
773 |
,NULL |
774 |
,NULL |
775 |
}; |
776 |
|
777 |
int ipopt_register(void){ |
778 |
return solver_register(&ipopt_internals); |
779 |
} |