Scippy

SCIP

Solving Constraint Integer Programs

nlpi_ipopt.cpp
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2020 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file nlpi_ipopt.cpp
17  * @ingroup NLPIS
18  * @brief Ipopt NLP interface
19  * @author Stefan Vigerske
20  * @author Benjamin Müller
21  *
22  * @todo warm starts
23  * @todo use new_x: Ipopt sets new_x = false if any function has been evaluated for the current x already, while oracle allows new_x to be false only if the current function has been evaluated for the current x before
24  * @todo influence output by SCIP verblevel, too, e.g., print strong warnings if SCIP verblevel is full; but currently we have no access to SCIP verblevel
25  * @todo if too few degrees of freedom, solve a slack-minimization problem instead?
26  *
27  * This file can only be compiled if Ipopt is available.
28  * Otherwise, to resolve public functions, use nlpi_ipopt_dummy.c.
29  * Since the dummy code is C instead of C++, it has been moved into a separate file.
30  */
31 
32 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
33 
34 #include "nlpi/nlpi_ipopt.h"
35 
36 #include "nlpi/nlpi.h"
37 #include "nlpi/nlpioracle.h"
38 #include "nlpi/exprinterpret.h"
39 #include "scip/interrupt.h"
40 #include "scip/pub_misc.h"
41 #include "scip/pub_message.h"
42 #include "scip/misc.h"
43 
44 #include <new> /* for std::bad_alloc */
45 #include <sstream>
46 
47 #ifdef __GNUC__
48 #pragma GCC diagnostic ignored "-Wshadow"
49 #endif
50 #include "IpoptConfig.h"
51 #include "IpIpoptApplication.hpp"
52 #include "IpIpoptCalculatedQuantities.hpp"
53 #include "IpSolveStatistics.hpp"
54 #include "IpJournalist.hpp"
55 #include "IpIpoptData.hpp"
56 #include "IpTNLPAdapter.hpp"
57 #include "IpOrigIpoptNLP.hpp"
58 #include "IpLapack.hpp"
59 #ifdef __GNUC__
60 #pragma GCC diagnostic warning "-Wshadow"
61 #endif
62 
63 #if (IPOPT_VERSION_MAJOR < 3 || (IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 12))
64 #error "The Ipopt interface requires at least 3.12."
65 #endif
66 
67 using namespace Ipopt;
68 
69 #define NLPI_NAME "ipopt" /**< short concise name of solver */
70 #define NLPI_DESC "Ipopt interface" /**< description of solver */
71 #define NLPI_PRIORITY 1000 /**< priority */
72 
73 #ifdef SCIP_DEBUG
74 #define DEFAULT_PRINTLEVEL J_WARNING /**< default print level of Ipopt */
75 #else
76 #define DEFAULT_PRINTLEVEL J_ERROR /**< default print level of Ipopt */
77 #endif
78 #define DEFAULT_MAXITER 3000 /**< default iteration limit for Ipopt */
79 
80 #define MAXPERTURB 0.01 /**< maximal perturbation of bounds in starting point heuristic */
81 #define FEASTOLFACTOR 0.05 /**< factor for user-given feasibility tolerance to get feasibility tolerance that is actually passed to Ipopt */
82 
83 #define DEFAULT_RANDSEED 71 /**< initial random seed */
84 
85 /* Convergence check (see ScipNLP::intermediate_callback)
86  *
87  * If the fastfail option is enabled, then we stop Ipopt if the reduction in
88  * primal infeasibility is not sufficient for a consecutive number of iterations.
89  * With the parameters as given below, we require Ipopt to
90  * - not increase the primal infeasibility after 5 iterations
91  * - reduce the primal infeasibility by at least 50% within 10 iterations
92  * - reduce the primal infeasibility by at least 90% within 30 iterations
93  * The targets are updated once they are reached and the limit on allowed iterations to reach the new target is reset.
94  *
95  * In certain situations, it is allowed to exceed an iteration limit:
96  * - If we are in the first 10 (convcheck_startiter) iterations.
97  * - If we are within 10 (convcheck_startiter) iterations after the restoration phase ended.
98  * The reason for this is that during feasibility restoration phase Ipopt aims completely on
99  * reducing constraint violation, completely forgetting the objective function.
100  * When returning from feasibility restoration and considering the original objective again,
101  * it is unlikely that Ipopt will continue to decrease primal infeasibility, since it may now target on
102  * more on optimality again. Thus, we do not check convergence for a number of iterations.
103  * - If the target on dual infeasibility reduction has been achieved, we are below twice the iteration limit, and
104  * we are not in restoration mode.
105  * The reason for this is that if Ipopt makes good progress towards optimality,
106  * we want to allow some more iterations where primal infeasibility is not reduced.
107  * However, in restoration mode, dual infeasibility does not correspond to the original problem and
108  * the complete aim is to restore primal infeasibility.
109  */
110 static const int convcheck_nchecks = 3; /**< number of convergence checks */
111 static const int convcheck_startiter = 10; /**< iteration where to start convergence checking */
112 static const int convcheck_maxiter[convcheck_nchecks] = { 5, 15, 30 }; /**< maximal number of iterations to achieve each convergence check */
113 static const SCIP_Real convcheck_minred[convcheck_nchecks] = { 1.0, 0.5, 0.1 }; /**< minimal required infeasibility reduction in each convergence check */
114 
115 class ScipNLP;
116 
117 struct SCIP_NlpiData
118 {
119 public:
120  BMS_BLKMEM* blkmem; /**< block memory */
121  SCIP_MESSAGEHDLR* messagehdlr; /**< message handler */
122  SCIP_Real infinity; /**< initial value for infinity */
123  std::string defoptions; /**< modified default options for Ipopt */
124 
125  /** constructor */
126  explicit SCIP_NlpiData(
127  BMS_BLKMEM* blkmem_ /**< block memory */
128  )
129  : blkmem(blkmem_),
130  messagehdlr(NULL),
132  { }
133 };
134 
135 struct SCIP_NlpiProblem
136 {
137 public:
138  SCIP_NLPIORACLE* oracle; /**< Oracle-helper to store and evaluate NLP */
139 
140  SmartPtr<IpoptApplication> ipopt; /**< Ipopt application */
141  SmartPtr<ScipNLP> nlp; /**< NLP in Ipopt form */
142  std::string optfile; /**< name of options file */
143  bool storeintermediate;/**< whether to store intermediate solutions */
144  bool fastfail; /**< whether to stop Ipopt if convergence seems slow */
145 
146  SCIP_Bool firstrun; /**< whether the next NLP solve will be the first one (with the current problem structure) */
147  SCIP_Real* initguess; /**< initial values for primal variables, or NULL if not known */
148 
149  SCIP_NLPSOLSTAT lastsolstat; /**< solution status from last run */
150  SCIP_NLPTERMSTAT lasttermstat; /**< termination status from last run */
151  SCIP_Real* lastsolprimals; /**< primal solution values from last run, if available */
152  SCIP_Real* lastsoldualcons; /**< dual solution values of constraints from last run, if available */
153  SCIP_Real* lastsoldualvarlb; /**< dual solution values of variable lower bounds from last run, if available */
154  SCIP_Real* lastsoldualvarub; /**< dual solution values of variable upper bounds from last run, if available */
155  SCIP_Real lastsolinfeas;/**< infeasibility (constraint violation) of solution stored in lastsolprimals */
156  int lastniter; /**< number of iterations in last run */
157  SCIP_Real lasttime; /**< time spend in last run */
158 
159  /** constructor */
161  : oracle(NULL),
162  storeintermediate(false), fastfail(false),
163  firstrun(TRUE), initguess(NULL),
164  lastsolstat(SCIP_NLPSOLSTAT_UNKNOWN), lasttermstat(SCIP_NLPTERMSTAT_OTHER),
165  lastsolprimals(NULL), lastsoldualcons(NULL), lastsoldualvarlb(NULL), lastsoldualvarub(NULL),
166  lastniter(-1), lasttime(-1.0)
167  { }
168 };
169 
170 /** TNLP implementation for SCIPs NLP */
171 class ScipNLP : public TNLP
172 {
173 private:
174  SCIP_NLPIPROBLEM* nlpiproblem; /**< NLPI problem data */
175  SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
176  BMS_BLKMEM* blkmem; /**< block memory */
177 
178  SCIP_Real conv_prtarget[convcheck_nchecks]; /**< target primal infeasibility for each convergence check */
179  SCIP_Real conv_dutarget[convcheck_nchecks]; /**< target dual infeasibility for each convergence check */
180  int conv_iterlim[convcheck_nchecks]; /**< iteration number where target primal infeasibility should to be achieved */
181  int conv_lastrestoiter; /**< last iteration number in restoration mode, or -1 if none */
182 
183 public:
184  bool approxhessian; /**< do we tell Ipopt to approximate the hessian? (may also be false if user set to approx. hessian via option file) */
185 
186  // cppcheck-suppress uninitMemberVar
187  /** constructor */
188  ScipNLP(
189  SCIP_NLPIPROBLEM* nlpiproblem_ = NULL,/**< NLPI problem data */
190  BMS_BLKMEM* blkmem_ = NULL /**< block memory */
191  )
192  : nlpiproblem(nlpiproblem_), randnumgen(NULL), blkmem(blkmem_), conv_lastrestoiter(-1), approxhessian(false)
193  {
194  assert(blkmem != NULL);
196  }
197 
198  /** destructor */
199  ~ScipNLP()
200  {
201  assert(randnumgen != NULL);
202  SCIPrandomFree(&randnumgen, blkmem);
203  }
204 
205  /** sets NLPI data structure */
206  void setNLPIPROBLEM(SCIP_NLPIPROBLEM* nlpiproblem_)
207  {
208  assert(nlpiproblem_ != NULL);
209  nlpiproblem = nlpiproblem_;
210  }
211 
212  /** Method to return some info about the nlp */
213  bool get_nlp_info(
214  Index& n, /**< place to store number of variables */
215  Index& m, /**< place to store number of constraints */
216  Index& nnz_jac_g, /**< place to store number of nonzeros in jacobian */
217  Index& nnz_h_lag, /**< place to store number of nonzeros in hessian */
218  IndexStyleEnum& index_style /**< place to store used index style (0-based or 1-based) */
219  );
220 
221  /** Method to return the bounds for my problem */
222  bool get_bounds_info(
223  Index n, /**< number of variables */
224  Number* x_l, /**< buffer to store lower bounds on variables */
225  Number* x_u, /**< buffer to store upper bounds on variables */
226  Index m, /**< number of constraints */
227  Number* g_l, /**< buffer to store lower bounds on constraints */
228  Number* g_u /**< buffer to store lower bounds on constraints */
229  );
230 
231  /** Method to return the starting point for the algorithm */
232  bool get_starting_point(
233  Index n, /**< number of variables */
234  bool init_x, /**< whether initial values for primal values are requested */
235  Number* x, /**< buffer to store initial primal values */
236  bool init_z, /**< whether initial values for dual values of variable bounds are requested */
237  Number* z_L, /**< buffer to store dual values for variable lower bounds */
238  Number* z_U, /**< buffer to store dual values for variable upper bounds */
239  Index m, /**< number of constraints */
240  bool init_lambda, /**< whether initial values for dual values of constraints are required */
241  Number* lambda /**< buffer to store dual values of constraints */
242  );
243 
244  /** Method to return the variables linearity. */
245  bool get_variables_linearity(
246  Index n, /**< number of variables */
247  LinearityType* var_types /**< buffer to store linearity types of variables */
248  );
249 
250  /** Method to return the constraint linearity. */
251  bool get_constraints_linearity(
252  Index m, /**< number of constraints */
253  LinearityType* const_types /**< buffer to store linearity types of constraints */
254  );
255 
256  /** Method to return the number of nonlinear variables. */
257  Index get_number_of_nonlinear_variables();
258 
259  /** Method to return the indices of the nonlinear variables */
260  bool get_list_of_nonlinear_variables(
261  Index num_nonlin_vars, /**< number of nonlinear variables */
262  Index* pos_nonlin_vars /**< array to fill with indices of nonlinear variables */
263  );
264 
265  /** Method to return metadata about variables and constraints */
266  bool get_var_con_metadata(
267  Index n, /**< number of variables */
268  StringMetaDataMapType& var_string_md, /**< variable meta data of string type */
269  IntegerMetaDataMapType& var_integer_md,/**< variable meta data of integer type */
270  NumericMetaDataMapType& var_numeric_md,/**< variable meta data of numeric type */
271  Index m, /**< number of constraints */
272  StringMetaDataMapType& con_string_md, /**< constraint meta data of string type */
273  IntegerMetaDataMapType& con_integer_md,/**< constraint meta data of integer type */
274  NumericMetaDataMapType& con_numeric_md /**< constraint meta data of numeric type */
275  );
276 
277  /** Method to return the objective value */
278  bool eval_f(
279  Index n, /**< number of variables */
280  const Number* x, /**< point to evaluate */
281  bool new_x, /**< whether some function evaluation method has been called for this point before */
282  Number& obj_value /**< place to store objective function value */
283  );
284 
285  /** Method to return the gradient of the objective */
286  bool eval_grad_f(
287  Index n, /**< number of variables */
288  const Number* x, /**< point to evaluate */
289  bool new_x, /**< whether some function evaluation method has been called for this point before */
290  Number* grad_f /**< buffer to store objective gradient */
291  );
292 
293  /** Method to return the constraint residuals */
294  bool eval_g(
295  Index n, /**< number of variables */
296  const Number* x, /**< point to evaluate */
297  bool new_x, /**< whether some function evaluation method has been called for this point before */
298  Index m, /**< number of constraints */
299  Number* g /**< buffer to store constraint function values */
300  );
301 
302  /** Method to return:
303  * 1) The structure of the jacobian (if "values" is NULL)
304  * 2) The values of the jacobian (if "values" is not NULL)
305  */
306  bool eval_jac_g(
307  Index n, /**< number of variables */
308  const Number* x, /**< point to evaluate */
309  bool new_x, /**< whether some function evaluation method has been called for this point before */
310  Index m, /**< number of constraints */
311  Index nele_jac, /**< number of nonzero entries in jacobian */
312  Index* iRow, /**< buffer to store row indices of nonzero jacobian entries, or NULL if values
313  * are requested */
314  Index* jCol, /**< buffer to store column indices of nonzero jacobian entries, or NULL if values
315  * are requested */
316  Number* values /**< buffer to store values of nonzero jacobian entries, or NULL if structure is
317  * requested */
318  );
319 
320  /** Method to return:
321  * 1) The structure of the hessian of the lagrangian (if "values" is NULL)
322  * 2) The values of the hessian of the lagrangian (if "values" is not NULL)
323  */
324  bool eval_h(
325  Index n, /**< number of variables */
326  const Number* x, /**< point to evaluate */
327  bool new_x, /**< whether some function evaluation method has been called for this point before */
328  Number obj_factor, /**< weight for objective function */
329  Index m, /**< number of constraints */
330  const Number* lambda, /**< weights for constraint functions */
331  bool new_lambda, /**< whether the hessian has been evaluated for these values of lambda before */
332  Index nele_hess, /**< number of nonzero entries in hessian */
333  Index* iRow, /**< buffer to store row indices of nonzero hessian entries, or NULL if values
334  * are requested */
335  Index* jCol, /**< buffer to store column indices of nonzero hessian entries, or NULL if values
336  * are requested */
337  Number* values /**< buffer to store values of nonzero hessian entries, or NULL if structure is requested */
338  );
339 
340  /** Method called by the solver at each iteration.
341  *
342  * Checks whether Ctrl-C was hit.
343  */
344  bool intermediate_callback(
345  AlgorithmMode mode, /**< current mode of algorithm */
346  Index iter, /**< current iteration number */
347  Number obj_value, /**< current objective value */
348  Number inf_pr, /**< current primal infeasibility */
349  Number inf_du, /**< current dual infeasibility */
350  Number mu, /**< current barrier parameter */
351  Number d_norm, /**< current gradient norm */
352  Number regularization_size,/**< current size of regularization */
353  Number alpha_du, /**< current dual alpha */
354  Number alpha_pr, /**< current primal alpha */
355  Index ls_trials, /**< current number of linesearch trials */
356  const IpoptData* ip_data, /**< pointer to Ipopt Data */
357  IpoptCalculatedQuantities* ip_cq /**< pointer to current calculated quantities */
358  );
359 
360  /** This method is called when the algorithm is complete so the TNLP can store/write the solution. */
361  void finalize_solution(
362  SolverReturn status, /**< solve and solution status */
363  Index n, /**< number of variables */
364  const Number* x, /**< primal solution values */
365  const Number* z_L, /**< dual values of variable lower bounds */
366  const Number* z_U, /**< dual values of variable upper bounds */
367  Index m, /**< number of constraints */
368  const Number* g, /**< values of constraints */
369  const Number* lambda, /**< dual values of constraints */
370  Number obj_value, /**< objective function value */
371  const IpoptData* data, /**< pointer to Ipopt Data */
372  IpoptCalculatedQuantities* cq /**< pointer to calculated quantities */
373  );
374 };
375 
376 /** A particular Ipopt::Journal implementation that uses the SCIP message routines for output. */
377 class ScipJournal : public Ipopt::Journal {
378 private:
379  /** reference to message handler pointer in NLPI data */
380  SCIP_MESSAGEHDLR*& messagehdlr;
381 
382 public:
383  ScipJournal(
384  const char* name, /**< name of journal */
385  Ipopt::EJournalLevel default_level, /**< default verbosity level */
386  SCIP_MESSAGEHDLR*& messagehdlr_ /**< pointer where to get message handler from */
387  )
388  : Ipopt::Journal(name, default_level),
389  messagehdlr(messagehdlr_)
390  { }
391 
392  ~ScipJournal() { }
393 
394 protected:
395  void PrintImpl(
396  Ipopt::EJournalCategory category, /**< category of message */
397  Ipopt::EJournalLevel level, /**< verbosity level of message */
398  const char* str /**< message to print */
399  )
400  {
401  if( level == J_ERROR )
402  {
404  }
405  else
406  {
407  SCIPmessagePrintInfo(messagehdlr, str);
408  }
409  }
410 
411  void PrintfImpl(
412  Ipopt::EJournalCategory category, /**< category of message */
413  Ipopt::EJournalLevel level, /**< verbosity level of message */
414  const char* pformat, /**< message printing format */
415  va_list ap /**< arguments of message */
416  )
417  {
418  if( level == J_ERROR )
419  {
420  SCIPmessageVPrintError(pformat, ap);
421  }
422  else
423  {
424  SCIPmessageVPrintInfo(messagehdlr, pformat, ap);
425  }
426  }
427 
428  void FlushBufferImpl() { }
429 };
430 
431 /** clears the last solution arrays and sets the solstat and termstat to unknown and other, resp. */
432 static
434  SCIP_NLPIPROBLEM* problem /**< data structure of problem */
435  )
436 {
437  assert(problem != NULL);
438 
445  problem->lastsolinfeas = SCIP_INVALID;
446 }
447 
448 /** sets feasibility tolerance parameters in Ipopt
449  *
450  * Sets tol and constr_viol_tol to FEASTOLFACTOR*feastol and acceptable_tol and acceptable_viol_tol to feastol.
451  * Since the users and Ipopts conception of feasibility may differ, we let Ipopt try to compute solutions
452  * that are more accurate (w.r.t. constraint violation) than requested by the user.
453  * Only if Ipopt has problems to achieve this accuracy, we also accept solutions that are accurate w.r.t. feastol only.
454  * The additional effort for computing a more accurate solution should be small if one can assume fast convergence when close to a local minimizer.
455  */
456 static
458  SCIP_NLPIPROBLEM* nlpiproblem,
459  SCIP_Real feastol
460  )
461 {
462  assert(nlpiproblem != NULL);
463 
464  nlpiproblem->ipopt->Options()->SetNumericValue("tol", FEASTOLFACTOR * feastol);
465  nlpiproblem->ipopt->Options()->SetNumericValue("constr_viol_tol", FEASTOLFACTOR * feastol);
466 
467  nlpiproblem->ipopt->Options()->SetNumericValue("acceptable_tol", feastol);
468  nlpiproblem->ipopt->Options()->SetNumericValue("acceptable_constr_viol_tol", feastol);
469 
470  /* It seem to be better to let Ipopt relax bounds a bit to ensure that a relative interior exists.
471  * However, if we relax the bounds too much, then the solutions tend to be slightly infeasible.
472  * If the user wants to set a tight feasibility tolerance, then (s)he has probably difficulties to compute accurate enough solutions.
473  * Thus, we turn off the bound_relax_factor completely if it would be below its default value of 1e-8.
474  */
475  nlpiproblem->ipopt->Options()->SetNumericValue("bound_relax_factor", feastol < 1e-8/FEASTOLFACTOR ? 0.0 : FEASTOLFACTOR * feastol);
476 }
477 
478 /** sets optimality tolerance parameters in Ipopt
479  *
480  * Sets dual_inf_tol and compl_inf_tol to opttol.
481  * We leave acceptable_dual_inf_tol and acceptable_compl_inf_tol untouched, which means that if Ipopt has convergence problems, then
482  * it can stop with a solution that is still feasible (see setFeastol), but essentially without a proof of local optimality.
483  * Note, that we report the solution as feasible only if Ipopt stopped on an "acceptable point" (see ScipNLP::finalize_solution).
484  *
485  * Note, that parameters tol and acceptable_tol (set in setFeastol) also apply.
486  */
487 static
489  SCIP_NLPIPROBLEM* nlpiproblem,
490  SCIP_Real opttol
491  )
492 {
493  assert(nlpiproblem != NULL);
494 
495  nlpiproblem->ipopt->Options()->SetNumericValue("dual_inf_tol", opttol);
496  nlpiproblem->ipopt->Options()->SetNumericValue("compl_inf_tol", opttol);
497 }
498 
499 /** copy method of NLP interface (called when SCIP copies plugins)
500  *
501  * input:
502  * - blkmem block memory of target SCIP
503  * - sourcenlpi the NLP interface to copy
504  * - targetnlpi buffer to store pointer to copy of NLP interface
505  */
506 static
507 SCIP_DECL_NLPICOPY(nlpiCopyIpopt)
508 {
509  SCIP_NLPIDATA* sourcedata;
510  SCIP_NLPIDATA* targetdata;
511 
512  assert(sourcenlpi != NULL);
513  assert(targetnlpi != NULL);
514 
515  sourcedata = SCIPnlpiGetData(sourcenlpi);
516  assert(sourcedata != NULL);
517 
518  SCIP_CALL( SCIPcreateNlpSolverIpopt(blkmem, targetnlpi) );
519  assert(*targetnlpi != NULL);
520 
521  SCIP_CALL( SCIPnlpiSetRealPar(*targetnlpi, NULL, SCIP_NLPPAR_INFINITY, sourcedata->infinity) );
522  SCIP_CALL( SCIPnlpiSetMessageHdlr(*targetnlpi, sourcedata->messagehdlr) );
523 
524  targetdata = SCIPnlpiGetData(*targetnlpi);
525  assert(targetdata != NULL);
526 
527  targetdata->defoptions = sourcedata->defoptions;
528 
529  return SCIP_OKAY;
530 }
531 
532 /** destructor of NLP interface to free nlpi data
533  *
534  * input:
535  * - nlpi datastructure for solver interface
536  */
537 static
538 SCIP_DECL_NLPIFREE(nlpiFreeIpopt)
539 {
540  SCIP_NLPIDATA* data;
541 
542  assert(nlpi != NULL);
543 
544  data = SCIPnlpiGetData(nlpi);
545  assert(data != NULL);
546 
547  delete data;
548 
549  return SCIP_OKAY;
550 }
551 
552 /** gets pointer for NLP solver to do dirty stuff
553  *
554  * input:
555  * - nlpi datastructure for solver interface
556  *
557  * return: void pointer to solver
558  */
559 static
560 SCIP_DECL_NLPIGETSOLVERPOINTER(nlpiGetSolverPointerIpopt)
561 {
562  assert(nlpi != NULL);
563 
564  return NULL;
565 }
566 
567 /** creates a problem instance
568  *
569  * input:
570  * - nlpi datastructure for solver interface
571  * - problem pointer to store the problem data
572  * - name name of problem, can be NULL
573  */
574 static
575 SCIP_DECL_NLPICREATEPROBLEM(nlpiCreateProblemIpopt)
576 {
577  SCIP_NLPIDATA* data;
578 
579  assert(nlpi != NULL);
580  assert(problem != NULL);
581 
582  data = SCIPnlpiGetData(nlpi);
583  assert(data != NULL);
584 
585  *problem = new SCIP_NLPIPROBLEM;
586  if( *problem == NULL )
587  return SCIP_NOMEMORY;
588 
589  SCIP_CALL( SCIPnlpiOracleCreate(data->blkmem, &(*problem)->oracle) );
590  SCIP_CALL( SCIPnlpiOracleSetInfinity((*problem)->oracle, data->infinity) );
591  SCIP_CALL( SCIPnlpiOracleSetProblemName((*problem)->oracle, name) );
592 
593  try
594  {
595  /* initialize IPOPT without default journal */
596  (*problem)->ipopt = new IpoptApplication(false);
597  if( IsNull((*problem)->ipopt) )
598  throw std::bad_alloc();
599 
600  /* plugin our journal to get output through SCIP message handler */
601  SmartPtr<Journal> jrnl = new ScipJournal("console", J_ITERSUMMARY, data->messagehdlr);
602  if( IsNull(jrnl) )
603  throw std::bad_alloc();
604  jrnl->SetPrintLevel(J_DBG, J_NONE);
605  if( !(*problem)->ipopt->Jnlst()->AddJournal(jrnl) )
606  {
607  SCIPerrorMessage("Failed to register ScipJournal for IPOPT output.");
608  }
609 
610  /* initialize Ipopt/SCIP NLP interface */
611  (*problem)->nlp = new ScipNLP(*problem, data->blkmem);
612  if( IsNull((*problem)->nlp) )
613  throw std::bad_alloc();
614  }
615  catch( const std::bad_alloc& )
616  {
617  SCIPerrorMessage("Not enough memory to initialize Ipopt.\n");
618  return SCIP_NOMEMORY;
619  }
620 
621  /* modify Ipopt's default settings to what we believe is appropriate */
622  (*problem)->ipopt->RegOptions()->AddStringOption2("store_intermediate", "whether to store the most feasible intermediate solutions", "no", "yes", "", "no", "", "useful when Ipopt looses a once found feasible solution and then terminates with an infeasible point");
623  (*problem)->ipopt->Options()->SetIntegerValue("print_level", DEFAULT_PRINTLEVEL);
624  /* (*problem)->ipopt->Options()->SetStringValue("print_timing_statistics", "yes"); */
625  (*problem)->ipopt->Options()->SetStringValue("mu_strategy", "adaptive");
626  (*problem)->ipopt->Options()->SetIntegerValue("max_iter", DEFAULT_MAXITER);
627  (*problem)->ipopt->Options()->SetNumericValue("nlp_lower_bound_inf", -data->infinity, false);
628  (*problem)->ipopt->Options()->SetNumericValue("nlp_upper_bound_inf", data->infinity, false);
629  (*problem)->ipopt->Options()->SetNumericValue("diverging_iterates_tol", data->infinity, false);
630  /* (*problem)->ipopt->Options()->SetStringValue("dependency_detector", "ma28"); */
631  /* if the expression interpreter does not give hessians, tell Ipopt to approximate hessian */
632  setFeastol(*problem, SCIP_DEFAULT_FEASTOL);
634 
635  /* apply user's given modifications to Ipopt's default settings */
636  if( data->defoptions.length() > 0 )
637  {
638  std::istringstream is(data->defoptions);
639 
640 #if (IPOPT_VERSION_MAJOR > 3) || (IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR > 12) || (IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR == 12 && IPOPT_VERSION_RELEASE >= 5)
641  if( !(*problem)->ipopt->Options()->ReadFromStream(*(*problem)->ipopt->Jnlst(), is, true) )
642 #else
643  if( !(*problem)->ipopt->Options()->ReadFromStream(*(*problem)->ipopt->Jnlst(), is) )
644 #endif
645  {
646  SCIPerrorMessage("Error when modifying Ipopt options using options string\n%s\n", data->defoptions.c_str());
647  return SCIP_ERROR;
648  }
649  }
650 
651  /* apply user's given options file (this one is NLPI problem specific) */
652  if( (*problem)->ipopt->Initialize((*problem)->optfile) != Solve_Succeeded )
653  {
654  SCIPerrorMessage("Error during initialization of Ipopt using optionfile \"%s\"\n", (*problem)->optfile.c_str());
655  return SCIP_ERROR;
656  }
657 
658 
659  return SCIP_OKAY;
660 }
661 
662 /** free a problem instance
663  *
664  * input:
665  * - nlpi datastructure for solver interface
666  * - problem pointer where problem data is stored
667  */
668 static
669 SCIP_DECL_NLPIFREEPROBLEM(nlpiFreeProblemIpopt)
670 {
671  assert(nlpi != NULL);
672  assert(problem != NULL);
673  assert(*problem != NULL);
674 
675  if( (*problem)->oracle != NULL )
676  {
677  SCIP_CALL( SCIPnlpiOracleFree(&(*problem)->oracle) );
678  }
679 
680  BMSfreeMemoryArrayNull(&(*problem)->initguess);
681  BMSfreeMemoryArrayNull(&(*problem)->lastsolprimals);
682  BMSfreeMemoryArrayNull(&(*problem)->lastsoldualcons);
683  BMSfreeMemoryArrayNull(&(*problem)->lastsoldualvarlb);
684  BMSfreeMemoryArrayNull(&(*problem)->lastsoldualvarub);
685 
686  delete *problem;
687  *problem = NULL;
688 
689  return SCIP_OKAY;
690 }
691 
692 /** gets pointer to solver-internal problem instance to do dirty stuff
693  *
694  * input:
695  * - nlpi datastructure for solver interface
696  * - problem datastructure for problem instance
697  *
698  * return: void pointer to problem instance
699  */
700 static
701 SCIP_DECL_NLPIGETPROBLEMPOINTER(nlpiGetProblemPointerIpopt)
702 {
703  assert(nlpi != NULL);
704  assert(problem != NULL);
705 
706  return GetRawPtr(problem->nlp);
707 }
708 
709 /** add variables
710  *
711  * input:
712  * - nlpi datastructure for solver interface
713  * - problem datastructure for problem instance
714  * - nvars number of variables
715  * - lbs lower bounds of variables, can be NULL if -infinity
716  * - ubs upper bounds of variables, can be NULL if +infinity
717  * - varnames names of variables, can be NULL
718  */
719 static
720 SCIP_DECL_NLPIADDVARS(nlpiAddVarsIpopt)
721 {
722  assert(nlpi != NULL);
723  assert(problem != NULL);
724  assert(problem->oracle != NULL);
725 
726  SCIP_CALL( SCIPnlpiOracleAddVars(problem->oracle, nvars, lbs, ubs, varnames) );
727 
728  problem->firstrun = TRUE;
729  BMSfreeMemoryArrayNull(&problem->initguess);
730  invalidateSolution(problem);
731 
732  return SCIP_OKAY;
733 }
734 
735 /** add constraints
736  *
737  * quadratic coefficiens: row oriented matrix for each constraint
738  *
739  * input:
740  * - nlpi datastructure for solver interface
741  * - problem datastructure for problem instance
742  * - ncons number of added constraints
743  * - lhss left hand sides of constraints
744  * - rhss right hand sides of constraints
745  * - linoffsets start index of each constraints linear coefficients in lininds and linvals
746  * length: ncons + 1, linoffsets[ncons] gives length of lininds and linvals
747  * may be NULL in case of no linear part
748  * - lininds variable indices
749  * may be NULL in case of no linear part
750  * - linvals coefficient values
751  * may be NULL in case of no linear part
752  * - nquadelems number of quadratic elements for each constraint
753  * may be NULL in case of no quadratic part
754  * - quadelems quadratic elements for each constraint
755  * may be NULL in case of no quadratic part
756  * - exprvaridxs indices of variables in expression tree, maps variable indices in expression tree to indices in nlp
757  * entry of array may be NULL in case of no expression tree
758  * may be NULL in case of no expression tree in any constraint
759  * - exprtrees expression tree for nonquadratic part of constraints
760  * entry of array may be NULL in case of no nonquadratic part
761  * may be NULL in case of no nonquadratic part in any constraint
762  * - names of constraints, may be NULL or entries may be NULL
763  */
764 static
765 SCIP_DECL_NLPIADDCONSTRAINTS(nlpiAddConstraintsIpopt)
766 {
767  assert(nlpi != NULL);
768  assert(problem != NULL);
769  assert(problem->oracle != NULL);
770 
771  SCIP_CALL( SCIPnlpiOracleAddConstraints(problem->oracle,
772  ncons, lhss, rhss,
773  nlininds, lininds, linvals,
774  nquadelems, quadelems,
775  exprvaridxs, exprtrees, names) );
776 
777  problem->firstrun = TRUE;
778  invalidateSolution(problem);
779 
780  return SCIP_OKAY;
781 }
782 
783 /** sets or overwrites objective, a minimization problem is expected
784  *
785  * May change sparsity pattern.
786  *
787  * input:
788  * - nlpi datastructure for solver interface
789  * - problem datastructure for problem instance
790  * - nlins number of linear variables
791  * - lininds variable indices
792  * may be NULL in case of no linear part
793  * - linvals coefficient values
794  * may be NULL in case of no linear part
795  * - nquadelems number of elements in matrix of quadratic part
796  * - quadelems elements of quadratic part
797  * may be NULL in case of no quadratic part
798  * - exprvaridxs indices of variables in expression tree, maps variable indices in expression tree to indices in nlp
799  * may be NULL in case of no expression tree
800  * - exprtree expression tree for nonquadratic part of objective function
801  * may be NULL in case of no nonquadratic part
802  * - constant objective value offset
803  */
804 static
805 SCIP_DECL_NLPISETOBJECTIVE(nlpiSetObjectiveIpopt)
806 {
807  assert(nlpi != NULL);
808  assert(problem != NULL);
809  assert(problem->oracle != NULL);
810 
811  /* We pass the objective gradient in dense form to Ipopt, so if the sparsity of that gradient changes, we do not need to reset Ipopt (firstrun=TRUE).
812  * However, if the sparsity of the Hessian matrix of the objective changes, then the sparsity pattern of the Hessian of the Lagrangian may change.
813  * Thus, reset Ipopt if the objective was and/or becomes nonlinear, but leave firstrun untouched if it was and stays linear.
814  */
815  if( nquadelems > 0 || exprtree != NULL || SCIPnlpiOracleGetConstraintDegree(problem->oracle, -1) > 1 )
816  problem->firstrun = TRUE;
817 
818  SCIP_CALL( SCIPnlpiOracleSetObjective(problem->oracle,
819  constant, nlins, lininds, linvals,
820  nquadelems, quadelems,
821  exprvaridxs, exprtree) );
822 
823  invalidateSolution(problem);
824 
825  return SCIP_OKAY;
826 }
827 
828 /** change variable bounds
829  *
830  * input:
831  * - nlpi datastructure for solver interface
832  * - problem datastructure for problem instance
833  * - nvars number of variables to change bounds
834  * - indices indices of variables to change bounds
835  * - lbs new lower bounds
836  * - ubs new upper bounds
837  */
838 static
839 SCIP_DECL_NLPICHGVARBOUNDS(nlpiChgVarBoundsIpopt)
840 {
841  int i;
842 
843  assert(nlpi != NULL);
844  assert(problem != NULL);
845  assert(problem->oracle != NULL);
846 
847  /* If some variable is fixed or unfixed, then better don't reoptimize the NLP in the next solve.
848  * Calling Optimize instead of ReOptimize should remove fixed variables from the problem that is solved by Ipopt.
849  * This way, the variable fixing is satisfied exactly in a solution, see also #1254.
850  */
851  for( i = 0; i < nvars && !problem->firstrun; ++i )
852  if( (SCIPnlpiOracleGetVarLbs(problem->oracle)[indices[i]] == SCIPnlpiOracleGetVarUbs(problem->oracle)[indices[i]]) != (lbs[i] == ubs[i]) )
853  problem->firstrun = TRUE;
854 
855  SCIP_CALL( SCIPnlpiOracleChgVarBounds(problem->oracle, nvars, indices, lbs, ubs) );
856 
857  invalidateSolution(problem);
858 
859  return SCIP_OKAY;
860 }
861 
862 /** change constraint bounds
863  *
864  * input:
865  * - nlpi datastructure for solver interface
866  * - problem datastructure for problem instance
867  * - nconss number of constraints to change sides
868  * - indices indices of constraints to change sides
869  * - lhss new left hand sides
870  * - rhss new right hand sides
871  */
872 static
873 SCIP_DECL_NLPICHGCONSSIDES(nlpiChgConsSidesIpopt)
874 {
875  assert(nlpi != NULL);
876  assert(problem != NULL);
877  assert(problem->oracle != NULL);
878 
879  SCIP_CALL( SCIPnlpiOracleChgConsSides(problem->oracle, nconss, indices, lhss, rhss) );
880 
881  invalidateSolution(problem);
882 
883  return SCIP_OKAY;
884 }
885 
886 /** delete a set of variables
887  *
888  * input:
889  * - nlpi datastructure for solver interface
890  * - problem datastructure for problem instance
891  * - nlpi datastructure for solver interface
892  * - dstats deletion status of vars; 1 if var should be deleted, 0 if not
893  *
894  * output:
895  * - dstats new position of var, -1 if var was deleted
896  */
897 static
898 SCIP_DECL_NLPIDELVARSET(nlpiDelVarSetIpopt)
899 {
900  assert(nlpi != NULL);
901  assert(problem != NULL);
902  assert(problem->oracle != NULL);
903 
904  SCIP_CALL( SCIPnlpiOracleDelVarSet(problem->oracle, dstats) );
905 
906  problem->firstrun = TRUE;
907  BMSfreeMemoryArrayNull(&problem->initguess); // @TODO keep initguess for remaining variables
908 
909  invalidateSolution(problem);
910 
911  return SCIP_OKAY;
912 }
913 
914 /** delete a set of constraints
915  *
916  * input:
917  * - nlpi datastructure for solver interface
918  * - problem datastructure for problem instance
919  * - dstats deletion status of rows; 1 if row should be deleted, 0 if not
920  *
921  * output:
922  * - dstats new position of row, -1 if row was deleted
923  */
924 static
925 SCIP_DECL_NLPIDELCONSSET(nlpiDelConstraintSetIpopt)
926 {
927  assert(nlpi != NULL);
928  assert(problem != NULL);
929  assert(problem->oracle != NULL);
930 
931  SCIP_CALL( SCIPnlpiOracleDelConsSet(problem->oracle, dstats) );
932 
933  problem->firstrun = TRUE;
934 
935  invalidateSolution(problem);
936 
937  return SCIP_OKAY;
938 }
939 
940 /** change one linear coefficient in a constraint or objective
941  *
942  * input:
943  * - nlpi datastructure for solver interface
944  * - problem datastructure for problem instance
945  * - idx index of constraint or -1 for objective
946  * - nvals number of values in linear constraint
947  * - varidxs indices of variable
948  * - vals new values for coefficient
949  *
950  * return: Error if coefficient did not exist before
951  */
952 static
953 SCIP_DECL_NLPICHGLINEARCOEFS(nlpiChgLinearCoefsIpopt)
954 {
955  assert(nlpi != NULL);
956  assert(problem != NULL);
957  assert(problem->oracle != NULL);
958 
959  SCIP_CALL( SCIPnlpiOracleChgLinearCoefs(problem->oracle, idx, nvals, varidxs, vals) );
960  invalidateSolution(problem);
961 
962  return SCIP_OKAY;
963 }
964 
965 /** change one coefficient in the quadratic part of a constraint or objective
966  *
967  * input:
968  * - nlpi datastructure for solver interface
969  * - problem datastructure for problem instance
970  * - idx index of constraint or -1 for objective
971  * - nquadelems number of entries in quadratic matrix to change
972  * - quadelems new elements in quadratic matrix (replacing already existing ones or adding new ones)
973  *
974  * return: Error if coefficient did not exist before
975  */
976 static
977 SCIP_DECL_NLPICHGQUADCOEFS(nlpiChgQuadraticCoefsIpopt)
978 {
979  assert(nlpi != NULL);
980  assert(problem != NULL);
981  assert(problem->oracle != NULL);
982 
983  SCIP_CALL( SCIPnlpiOracleChgQuadCoefs(problem->oracle, idx, nquadelems, quadelems) );
984  invalidateSolution(problem);
985 
986  return SCIP_OKAY;
987 }
988 
989 /** replaces the expression tree of a constraint or objective
990  *
991  * input:
992  * - nlpi datastructure for solver interface
993  * - problem datastructure for problem instance
994  * - idxcons index of constraint or -1 for objective
995  * - exprtree new expression tree for constraint or objective, or NULL to only remove previous tree
996  */
997 static
998 SCIP_DECL_NLPICHGEXPRTREE(nlpiChgExprtreeIpopt)
999 {
1000  assert(nlpi != NULL);
1001  assert(problem != NULL);
1002  assert(problem->oracle != NULL);
1003 
1004  SCIP_CALL( SCIPnlpiOracleChgExprtree(problem->oracle, idxcons, exprvaridxs, exprtree) );
1005  invalidateSolution(problem);
1006 
1007  return SCIP_OKAY;
1008 }
1009 
1010 /** change the value of one parameter in the nonlinear part
1011  *
1012  * input:
1013  * - nlpi datastructure for solver interface
1014  * - problem datastructure for problem instance
1015  * - idxcons index of constraint or -1 for objective
1016  * - idxparam index of parameter
1017  * - value new value for nonlinear parameter
1018  *
1019  * return: Error if parameter does not exist
1020  */
1021 static
1022 SCIP_DECL_NLPICHGNONLINCOEF(nlpiChgNonlinCoefIpopt)
1023 {
1024  assert(nlpi != NULL);
1025  assert(problem != NULL);
1026  assert(problem->oracle != NULL);
1027 
1028  SCIP_CALL( SCIPnlpiOracleChgExprParam(problem->oracle, idxcons, idxparam, value) );
1029  invalidateSolution(problem);
1030 
1031  return SCIP_OKAY;
1032 }
1033 
1034 /** change the constant offset in the objective
1035  *
1036  * input:
1037  * - nlpi datastructure for solver interface
1038  * - problem datastructure for problem instance
1039  * - objconstant new value for objective constant
1040  */
1041 static
1042 SCIP_DECL_NLPICHGOBJCONSTANT( nlpiChgObjConstantIpopt )
1043 {
1044  assert(nlpi != NULL);
1045  assert(problem != NULL);
1046  assert(problem->oracle != NULL);
1047 
1048  SCIP_CALL( SCIPnlpiOracleChgObjConstant(problem->oracle, objconstant) );
1049 
1050  return SCIP_OKAY;
1051 }
1052 
1053 /** sets initial guess for primal variables
1054  *
1055  * input:
1056  * - nlpi datastructure for solver interface
1057  * - problem datastructure for problem instance
1058  * - primalvalues initial primal values for variables, or NULL to clear previous values
1059  * - consdualvalues initial dual values for constraints, or NULL to clear previous values
1060  * - varlbdualvalues initial dual values for variable lower bounds, or NULL to clear previous values
1061  * - varubdualvalues initial dual values for variable upper bounds, or NULL to clear previous values
1062  */
1063 static
1064 SCIP_DECL_NLPISETINITIALGUESS(nlpiSetInitialGuessIpopt)
1065 {
1066  assert(nlpi != NULL);
1067  assert(problem != NULL);
1068  assert(problem->oracle != NULL);
1069 
1070  if( primalvalues != NULL )
1071  {
1072  if( !problem->initguess )
1073  {
1074  if( BMSduplicateMemoryArray(&problem->initguess, primalvalues, SCIPnlpiOracleGetNVars(problem->oracle)) == NULL )
1075  return SCIP_NOMEMORY;
1076  }
1077  else
1078  {
1079  BMScopyMemoryArray(problem->initguess, primalvalues, SCIPnlpiOracleGetNVars(problem->oracle));
1080  }
1081  }
1082  else
1083  {
1084  BMSfreeMemoryArrayNull(&problem->initguess);
1085  }
1086 
1087  return SCIP_OKAY;
1088 }
1089 
1090 /** tries to solve NLP
1091  *
1092  * input:
1093  * - nlpi datastructure for solver interface
1094  * - problem datastructure for problem instance
1095  */
1096 static
1097 SCIP_DECL_NLPISOLVE(nlpiSolveIpopt)
1098 {
1099  ApplicationReturnStatus status;
1100 
1101  assert(nlpi != NULL);
1102  assert(problem != NULL);
1103  assert(problem->oracle != NULL);
1104 
1105  assert(IsValid(problem->ipopt));
1106  assert(IsValid(problem->nlp));
1107 
1108  problem->nlp->setNLPIPROBLEM(problem);
1109 
1110  problem->lastniter = -1;
1111  problem->lasttime = -1.0;
1112  problem->lastsolinfeas = SCIP_INVALID;
1113 
1114  try
1115  {
1116  SmartPtr<SolveStatistics> stats;
1117 
1118  if( problem->firstrun )
1119  {
1121 
1122  cap = SCIPexprintGetCapability() & SCIPnlpiOracleGetEvalCapability(problem->oracle);
1123 
1124  /* if the expression interpreter or some user expression do not support function values and gradients and Hessians, and the problem is not at most quadratic,
1125  * change NLP parameters or give an error
1126  */
1128  && SCIPnlpiOracleGetMaxDegree(problem->oracle) > 2 )
1129  {
1130  /* @todo could enable Jacobian approximation in Ipopt */
1133  {
1134  SCIPerrorMessage("Do not have expression interpreter that can compute function values and gradients. Cannot solve NLP with Ipopt.\n");
1135  problem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
1136  problem->lasttermstat = SCIP_NLPTERMSTAT_OTHER;
1137  return SCIP_OKAY;
1138  }
1139 
1140  /* enable Hessian approximation if we are nonquadratic and the expression interpreter or user expression do not support Hessians */
1141  if( !(cap & SCIP_EXPRINTCAPABILITY_HESSIAN) )
1142  {
1143  problem->ipopt->Options()->SetStringValue("hessian_approximation", "limited-memory");
1144  problem->nlp->approxhessian = true;
1145  }
1146  else
1147  problem->nlp->approxhessian = false;
1148  }
1149 
1150 #ifdef SCIP_DEBUG
1151  problem->ipopt->Options()->SetStringValue("derivative_test", problem->nlp->approxhessian ? "first-order" : "second-order");
1152 #endif
1153 
1154  status = problem->ipopt->OptimizeTNLP(GetRawPtr(problem->nlp));
1155  }
1156  else
1157  {
1158  status = problem->ipopt->ReOptimizeTNLP(GetRawPtr(problem->nlp));
1159  }
1160 
1161  // catch the very bad status codes
1162  switch( status ) {
1163  case Invalid_Problem_Definition:
1164  case Invalid_Option:
1165  case Unrecoverable_Exception:
1166  case NonIpopt_Exception_Thrown:
1167  SCIPerrorMessage("Ipopt returned with application return status %d\n", status);
1168  return SCIP_ERROR;
1169  case Internal_Error:
1170  // could be a fail in the linear solver
1171  SCIPerrorMessage("Ipopt returned with status \"Internal Error\"\n");
1172  invalidateSolution(problem);
1173  problem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
1174  problem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
1175  break;
1176  case Insufficient_Memory:
1177  SCIPerrorMessage("Ipopt returned with status \"Insufficient Memory\"\n");
1178  return SCIP_NOMEMORY;
1179  case Invalid_Number_Detected:
1180  SCIPdebugMessage("Ipopt failed because of an invalid number in function or derivative value\n");
1181  invalidateSolution(problem);
1182  problem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
1183  problem->lasttermstat = SCIP_NLPTERMSTAT_EVALERR;
1184  break;
1185  default: ;
1186  }
1187 
1188  stats = problem->ipopt->Statistics();
1189  if( IsValid(stats) )
1190  {
1191  problem->lastniter = stats->IterationCount();
1192  problem->lasttime = stats->TotalCPUTime();
1193  }
1194  else
1195  {
1196  /* Ipopt does not provide access to the statistics when all variables have been fixed */
1197  problem->lastniter = 0;
1198  problem->lasttime = 0.0;
1199  }
1200  }
1201  catch( IpoptException& except )
1202  {
1203  SCIPerrorMessage("Ipopt returned with exception: %s\n", except.Message().c_str());
1204  return SCIP_ERROR;
1205  }
1206 
1207  problem->firstrun = FALSE;
1208 
1209  return SCIP_OKAY;
1210 }
1211 
1212 /** gives solution status
1213  *
1214  * input:
1215  * - nlpi datastructure for solver interface
1216  * - problem datastructure for problem instance
1217  *
1218  * return: Solution Status
1219  */
1220 static
1221 SCIP_DECL_NLPIGETSOLSTAT(nlpiGetSolstatIpopt)
1222 {
1223  assert(nlpi != NULL);
1224  assert(problem != NULL);
1225 
1226  return problem->lastsolstat;
1227 }
1228 
1229 /** gives termination reason
1230  *
1231  * input:
1232  * - nlpi datastructure for solver interface
1233  * - problem datastructure for problem instance
1234  *
1235  * return: Termination Status
1236  */
1237 static
1238 SCIP_DECL_NLPIGETTERMSTAT(nlpiGetTermstatIpopt)
1239 {
1240  assert(nlpi != NULL);
1241  assert(problem != NULL);
1242 
1243  return problem->lasttermstat;
1244 }
1245 
1246 /** gives primal and dual solution values
1247  *
1248  * input:
1249  * - nlpi datastructure for solver interface
1250  * - problem datastructure for problem instance
1251  * - primalvalues buffer to store pointer to array to primal values, or NULL if not needed
1252  * - consdualvalues buffer to store pointer to array to dual values of constraints, or NULL if not needed
1253  * - varlbdualvalues buffer to store pointer to array to dual values of variable lower bounds, or NULL if not needed
1254  * - varubdualvalues buffer to store pointer to array to dual values of variable lower bounds, or NULL if not needed
1255  * - objval buffer store the objective value, or NULL if not needed
1256  */
1257 static
1258 SCIP_DECL_NLPIGETSOLUTION(nlpiGetSolutionIpopt)
1259 {
1260  assert(nlpi != NULL);
1261  assert(problem != NULL);
1262 
1263  if( primalvalues != NULL )
1264  *primalvalues = problem->lastsolprimals;
1265 
1266  if( consdualvalues != NULL )
1267  *consdualvalues = problem->lastsoldualcons;
1268 
1269  if( varlbdualvalues != NULL )
1270  *varlbdualvalues = problem->lastsoldualvarlb;
1271 
1272  if( varubdualvalues != NULL )
1273  *varubdualvalues = problem->lastsoldualvarub;
1274 
1275  if( objval != NULL )
1276  {
1277  if( problem->lastsolprimals != NULL )
1278  {
1279  /* TODO store last solution value instead of reevaluating the objective function */
1280  SCIP_CALL( SCIPnlpiOracleEvalObjectiveValue(problem->oracle, problem->lastsolprimals, objval) );
1281  }
1282  else
1283  *objval = SCIP_INVALID;
1284  }
1285 
1286  return SCIP_OKAY;
1287 }
1288 
1289 /** gives solve statistics
1290  *
1291  * input:
1292  * - nlpi datastructure for solver interface
1293  * - problem datastructure for problem instance
1294  * - statistics pointer to store statistics
1295  *
1296  * output:
1297  * - statistics solve statistics
1298  */
1299 static
1300 SCIP_DECL_NLPIGETSTATISTICS(nlpiGetStatisticsIpopt)
1301 {
1302  assert(nlpi != NULL);
1303  assert(problem != NULL);
1304 
1305  SCIPnlpStatisticsSetNIterations(statistics, problem->lastniter);
1306  SCIPnlpStatisticsSetTotalTime (statistics, problem->lasttime);
1307 
1308  return SCIP_OKAY;
1309 }
1310 
1311 /** gives required size of a buffer to store a warmstart object
1312  *
1313  * input:
1314  * - nlpi datastructure for solver interface
1315  * - problem datastructure for problem instance
1316  * - size pointer to store required size for warmstart buffer
1317  *
1318  * output:
1319  * - size required size for warmstart buffer
1320  */
1321 static
1322 SCIP_DECL_NLPIGETWARMSTARTSIZE(nlpiGetWarmstartSizeIpopt)
1323 {
1324  SCIPerrorMessage("method of Ipopt nonlinear solver is not implemented\n");
1325  return SCIP_ERROR;
1326 }
1327 
1328 /** stores warmstart information in buffer
1329  *
1330  * Required size of buffer should have been obtained by SCIPnlpiGetWarmstartSize before.
1331  *
1332  * input:
1333  * - nlpi datastructure for solver interface
1334  * - problem datastructure for problem instance
1335  * - buffer memory to store warmstart information
1336  *
1337  * output:
1338  * - buffer warmstart information in solver specific data structure
1339  */
1340 static
1341 SCIP_DECL_NLPIGETWARMSTARTMEMO(nlpiGetWarmstartMemoIpopt)
1342 {
1343  SCIPerrorMessage("method of Ipopt nonlinear solver is not implemented\n");
1344  return SCIP_ERROR;
1345 }
1346 
1347 /** sets warmstart information in solver
1348  *
1349  * Write warmstart to buffer.
1350  *
1351  * input:
1352  * - nlpi datastructure for solver interface
1353  * - problem datastructure for problem instance
1354  * - buffer warmstart information
1355  */
1356 static
1357 SCIP_DECL_NLPISETWARMSTARTMEMO(nlpiSetWarmstartMemoIpopt)
1358 {
1359  SCIPerrorMessage("method of Ipopt nonlinear solver is not implemented\n");
1360  SCIPABORT();
1361  return SCIP_OKAY;
1362 }
1363 
1364 /** gets integer parameter of NLP
1365  *
1366  * input:
1367  * - nlpi datastructure for solver interface
1368  * - problem datastructure for problem instance
1369  * - type parameter number
1370  * - ival pointer to store the parameter value
1371  *
1372  * output:
1373  * - ival parameter value
1374  */
1375 static
1376 SCIP_DECL_NLPIGETINTPAR(nlpiGetIntParIpopt)
1377 {
1378  assert(nlpi != NULL);
1379  assert(ival != NULL);
1380  assert(problem != NULL);
1381  assert(IsValid(problem->ipopt));
1382 
1383  //@TODO try-catch block for Ipopt exceptions
1384  switch( type )
1385  {
1387  {
1388  *ival = 1;
1389  break;
1390  }
1391 
1392  case SCIP_NLPPAR_VERBLEVEL:
1393  {
1394  int printlevel;
1395  problem->ipopt->Options()->GetIntegerValue("print_level", printlevel, "");
1396  if( printlevel <= J_ERROR )
1397  *ival = 0;
1398  else if( printlevel >= J_DETAILED )
1399  *ival = printlevel - J_ITERSUMMARY + 1;
1400  else /* J_SUMMARY or J_WARNING or J_ITERSUMMARY */
1401  *ival = 1;
1402  break;
1403  }
1404 
1405  case SCIP_NLPPAR_FEASTOL:
1406  {
1407  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
1408  return SCIP_PARAMETERWRONGTYPE;
1409  }
1410 
1411  case SCIP_NLPPAR_RELOBJTOL:
1412  {
1413  SCIPerrorMessage("relative objective tolerance parameter is of type real.\n");
1414  return SCIP_PARAMETERWRONGTYPE;
1415  }
1416 
1417  case SCIP_NLPPAR_LOBJLIM:
1418  {
1419  SCIPerrorMessage("objective limit parameter is of type real.\n");
1420  return SCIP_PARAMETERWRONGTYPE;
1421  }
1422 
1423  case SCIP_NLPPAR_INFINITY:
1424  {
1425  SCIPerrorMessage("infinity parameter is of type real.\n");
1426  return SCIP_PARAMETERWRONGTYPE;
1427  }
1428 
1429  case SCIP_NLPPAR_ITLIM:
1430  {
1431  problem->ipopt->Options()->GetIntegerValue("max_iter", *ival, "");
1432  break;
1433  }
1434 
1435  case SCIP_NLPPAR_TILIM:
1436  {
1437  SCIPerrorMessage("time limit parameter is of type real.\n");
1438  return SCIP_PARAMETERWRONGTYPE;
1439  }
1440 
1441  case SCIP_NLPPAR_OPTFILE:
1442  {
1443  SCIPerrorMessage("optfile parameter is of type string.\n");
1444  return SCIP_PARAMETERWRONGTYPE;
1445  }
1446 
1447  case SCIP_NLPPAR_FASTFAIL:
1448  {
1449  *ival = problem->fastfail ? 1 : 0;
1450  break;
1451  }
1452 
1453  default:
1454  {
1455  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
1456  return SCIP_PARAMETERUNKNOWN;
1457  }
1458  }
1459 
1460  return SCIP_OKAY;
1461 }
1462 
1463 /** sets integer parameter of NLP
1464  *
1465  * input:
1466  * - nlpi datastructure for solver interface
1467  * - problem datastructure for problem instance
1468  * - type parameter number
1469  * - ival parameter value
1470  */
1471 static
1472 SCIP_DECL_NLPISETINTPAR(nlpiSetIntParIpopt)
1473 {
1474  assert(nlpi != NULL);
1475  assert(problem != NULL);
1476  assert(IsValid(problem->ipopt));
1477 
1478  switch( type )
1479  {
1481  {
1482  if( ival == 0 || ival == 1 )
1483  {
1484  SCIP_NLPIDATA* data;
1485 
1486  data = SCIPnlpiGetData(nlpi);
1487  assert(data != NULL);
1488 
1489  SCIPmessagePrintWarning(data->messagehdlr, "from scratch parameter not supported by Ipopt interface yet. Ignored.\n");
1490  }
1491  else
1492  {
1493  SCIPerrorMessage("Value %d for parameter from scratch out of range {0, 1}\n", ival);
1494  return SCIP_PARAMETERWRONGVAL;
1495  }
1496  break;
1497  }
1498 
1499  case SCIP_NLPPAR_VERBLEVEL:
1500  {
1501  switch( ival )
1502  {
1503  case 0:
1504  problem->ipopt->Options()->SetIntegerValue("print_level", J_ERROR);
1505  break;
1506  case 1:
1507  problem->ipopt->Options()->SetIntegerValue("print_level", J_ITERSUMMARY);
1508  break;
1509  case 2:
1510  problem->ipopt->Options()->SetIntegerValue("print_level", J_DETAILED);
1511  break;
1512  default:
1513  if( ival > 2 )
1514  {
1515  problem->ipopt->Options()->SetIntegerValue("print_level", MIN(J_ITERSUMMARY + (ival-1), J_ALL));
1516  break;
1517  }
1518  else
1519  {
1520  SCIPerrorMessage("Value %d for parameter from verbosity level out of range {0, 1, 2}\n", ival);
1521  return SCIP_PARAMETERWRONGVAL;
1522  }
1523  }
1524  break;
1525  }
1526 
1527  case SCIP_NLPPAR_FEASTOL:
1528  {
1529  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
1530  return SCIP_PARAMETERWRONGTYPE;
1531  }
1532 
1533  case SCIP_NLPPAR_RELOBJTOL:
1534  {
1535  SCIPerrorMessage("relative objective tolerance parameter is of type real.\n");
1536  return SCIP_PARAMETERWRONGTYPE;
1537  }
1538 
1539  case SCIP_NLPPAR_LOBJLIM:
1540  {
1541  SCIPerrorMessage("objective limit parameter is of type real.\n");
1542  return SCIP_PARAMETERWRONGTYPE;
1543  }
1544 
1545  case SCIP_NLPPAR_INFINITY:
1546  {
1547  SCIPerrorMessage("infinity parameter is of type real.\n");
1548  return SCIP_PARAMETERWRONGTYPE;
1549  }
1550 
1551  case SCIP_NLPPAR_ITLIM:
1552  {
1553  if( ival >= 0 )
1554  {
1555  problem->ipopt->Options()->SetIntegerValue("max_iter", ival);
1556  }
1557  else
1558  {
1559  SCIPerrorMessage("Value %d for parameter iteration limit is negative\n", ival);
1560  return SCIP_PARAMETERWRONGVAL;
1561  }
1562  break;
1563  }
1564 
1565  case SCIP_NLPPAR_TILIM:
1566  {
1567  SCIPerrorMessage("time limit parameter is of type real.\n");
1568  return SCIP_PARAMETERWRONGTYPE;
1569  }
1570 
1571  case SCIP_NLPPAR_OPTFILE:
1572  {
1573  SCIPerrorMessage("optfile parameter is of type string.\n");
1574  return SCIP_PARAMETERWRONGTYPE;
1575  }
1576 
1577  case SCIP_NLPPAR_FASTFAIL:
1578  {
1579  if( ival == 0 || ival == 1 )
1580  {
1581  problem->fastfail = (bool)ival;
1582  problem->storeintermediate = (bool)ival;
1583  }
1584  else
1585  {
1586  SCIPerrorMessage("Value %d for parameter fastfail out of range {0, 1}\n", ival);
1587  return SCIP_PARAMETERWRONGVAL;
1588  }
1589  break;
1590  }
1591 
1592  default:
1593  {
1594  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
1595  return SCIP_PARAMETERUNKNOWN;
1596  }
1597  }
1598 
1599  return SCIP_OKAY;
1600 }
1601 
1602 /** gets floating point parameter of NLP
1603  *
1604  * input:
1605  * - nlpi datastructure for solver interface
1606  * - problem datastructure for problem instance, can be NULL only if type == SCIP_NLPPAR_INFINITY
1607  * - type parameter number
1608  * - dval pointer to store the parameter value
1609  *
1610  * output:
1611  * - dval parameter value
1612  */
1613 static
1614 SCIP_DECL_NLPIGETREALPAR(nlpiGetRealParIpopt)
1615 {
1616  assert(nlpi != NULL);
1617  assert(dval != NULL);
1618  assert(type == SCIP_NLPPAR_INFINITY || problem != NULL);
1619  assert(type == SCIP_NLPPAR_INFINITY || IsValid(problem->ipopt));
1620 
1621  switch( type )
1622  {
1624  {
1625  SCIPerrorMessage("from scratch parameter is of type int.\n");
1626  return SCIP_PARAMETERWRONGTYPE;
1627  }
1628 
1629  case SCIP_NLPPAR_VERBLEVEL:
1630  {
1631  SCIPerrorMessage("verbosity level parameter is of type int.\n");
1632  return SCIP_PARAMETERWRONGTYPE;
1633  }
1634 
1635  case SCIP_NLPPAR_FEASTOL:
1636  {
1637  problem->ipopt->Options()->GetNumericValue("acceptable_constr_viol_tol", *dval, "");
1638  break;
1639  }
1640 
1641  case SCIP_NLPPAR_RELOBJTOL:
1642  {
1643  problem->ipopt->Options()->GetNumericValue("dual_inf_tol", *dval, "");
1644  break;
1645  }
1646 
1647  case SCIP_NLPPAR_LOBJLIM:
1648  {
1649  *dval = -SCIPnlpiOracleGetInfinity(problem->oracle);
1650  break;
1651  }
1652 
1653  case SCIP_NLPPAR_INFINITY:
1654  {
1655  if( problem )
1656  {
1657  *dval = SCIPnlpiOracleGetInfinity(problem->oracle);
1658  }
1659  else
1660  {
1661  SCIP_NLPIDATA* data = SCIPnlpiGetData(nlpi);
1662  assert(data != NULL);
1663  *dval = data->infinity;
1664  }
1665  break;
1666  }
1667 
1668  case SCIP_NLPPAR_ITLIM:
1669  {
1670  SCIPerrorMessage("iteration limit parameter is of type int.\n");
1671  return SCIP_PARAMETERWRONGTYPE;
1672  }
1673 
1674  case SCIP_NLPPAR_TILIM:
1675  {
1676  problem->ipopt->Options()->GetNumericValue("max_cpu_time", *dval, "");
1677  break;
1678  }
1679 
1680  case SCIP_NLPPAR_OPTFILE:
1681  {
1682  SCIPerrorMessage("option file parameter is of type string.\n");
1683  return SCIP_PARAMETERWRONGTYPE;
1684  }
1685 
1686  case SCIP_NLPPAR_FASTFAIL:
1687  {
1688  SCIPerrorMessage("fastfail parameter is of type int.\n");
1689  return SCIP_PARAMETERWRONGTYPE;
1690  }
1691 
1692  default:
1693  {
1694  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
1695  return SCIP_PARAMETERUNKNOWN;
1696  }
1697  }
1698 
1699  return SCIP_OKAY;
1700 }
1701 
1702 /** sets floating point parameter of NLP
1703  *
1704  * input:
1705  * - nlpi datastructure for solver interface
1706  * - problem datastructure for problem instance, can be NULL only if type == SCIP_NLPPAR_INFINITY
1707  * - type parameter number
1708  * - dval parameter value
1709  */
1710 static
1711 SCIP_DECL_NLPISETREALPAR(nlpiSetRealParIpopt)
1712 {
1713  assert(nlpi != NULL);
1714  assert(type == SCIP_NLPPAR_INFINITY || problem != NULL);
1715  assert(type == SCIP_NLPPAR_INFINITY || IsValid(problem->ipopt));
1716 
1717  switch( type )
1718  {
1720  {
1721  SCIPerrorMessage("from scratch parameter is of type int.\n");
1722  return SCIP_PARAMETERWRONGTYPE;
1723  }
1724 
1725  case SCIP_NLPPAR_VERBLEVEL:
1726  {
1727  SCIPerrorMessage("verbosity level parameter is of type int.\n");
1728  return SCIP_PARAMETERWRONGTYPE;
1729  }
1730 
1731  case SCIP_NLPPAR_FEASTOL:
1732  {
1733  if( dval >= 0 )
1734  {
1735  setFeastol(problem, dval);
1736  }
1737  else
1738  {
1739  SCIPerrorMessage("Value %g for parameter feasibility tolerance is negative\n", dval);
1740  return SCIP_PARAMETERWRONGVAL;
1741  }
1742  break;
1743  }
1744 
1745  case SCIP_NLPPAR_RELOBJTOL:
1746  {
1747  if( dval >= 0 )
1748  {
1749  setOpttol(problem, dval);
1750  }
1751  else
1752  {
1753  SCIPerrorMessage("Value %g for parameter relative objective tolerance is negative\n", dval);
1754  return SCIP_PARAMETERWRONGVAL;
1755  }
1756  break;
1757  }
1758 
1759  case SCIP_NLPPAR_LOBJLIM:
1760  {
1761  SCIP_NLPIDATA* data;
1762 
1763  data = SCIPnlpiGetData(nlpi);
1764  assert(data != NULL);
1765 
1766  SCIPmessagePrintWarning(data->messagehdlr, "Parameter lower objective limit not supported by Ipopt interface yet. Ignored.\n");
1767  break;
1768  }
1769 
1770  case SCIP_NLPPAR_INFINITY:
1771  {
1772  if( dval < 0.0 )
1773  return SCIP_PARAMETERWRONGVAL;
1774  if( problem )
1775  {
1776  problem->ipopt->Options()->SetNumericValue("diverging_iterates_tol", dval);
1777  problem->ipopt->Options()->SetNumericValue("nlp_lower_bound_inf", -dval);
1778  problem->ipopt->Options()->SetNumericValue("nlp_upper_bound_inf", dval);
1779  SCIPnlpiOracleSetInfinity(problem->oracle, dval);
1780  }
1781  else
1782  {
1783  SCIP_NLPIDATA* data = SCIPnlpiGetData(nlpi);
1784  assert(data != NULL);
1785  data->infinity = dval;
1786  }
1787  break;
1788  }
1789 
1790  case SCIP_NLPPAR_ITLIM:
1791  {
1792  SCIPerrorMessage("iteration limit parameter is of type int.\n");
1793  return SCIP_PARAMETERWRONGTYPE;
1794  }
1795 
1796  case SCIP_NLPPAR_TILIM:
1797  {
1798  if( dval >= 0.0 )
1799  {
1800  /* Ipopt doesn't like a setting of exactly 0 for the max_cpu_time, so increase as little as possible in that case */
1801  problem->ipopt->Options()->SetNumericValue("max_cpu_time", MAX(dval, DBL_MIN));
1802  }
1803  else
1804  {
1805  SCIPerrorMessage("Value %g for parameter time limit is negative\n", dval);
1806  return SCIP_PARAMETERWRONGVAL;
1807  }
1808  break;
1809  }
1810 
1811  case SCIP_NLPPAR_OPTFILE:
1812  {
1813  SCIPerrorMessage("option file parameter is of type string.\n");
1814  return SCIP_PARAMETERWRONGTYPE;
1815  }
1816 
1817  case SCIP_NLPPAR_FASTFAIL:
1818  {
1819  SCIPerrorMessage("fastfail parameter is of type int.\n");
1820  return SCIP_PARAMETERWRONGTYPE;
1821  }
1822 
1823  default:
1824  {
1825  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
1826  return SCIP_PARAMETERUNKNOWN;
1827  }
1828  }
1829 
1830  return SCIP_OKAY;
1831 }
1832 
1833 /** gets string parameter of NLP
1834  *
1835  * input:
1836  * - nlpi NLP interface structure
1837  * - problem datastructure for problem instance
1838  * - type parameter number
1839  * - sval pointer to store the string value, the user must not modify the string
1840  *
1841  * output:
1842  * - sval parameter value
1843  */
1844 static
1845 SCIP_DECL_NLPIGETSTRINGPAR( nlpiGetStringParIpopt )
1846 {
1847  assert(nlpi != NULL);
1848  assert(problem != NULL);
1849 
1850  switch( type )
1851  {
1853  {
1854  SCIPerrorMessage("from scratch parameter is of type int.\n");
1855  return SCIP_PARAMETERWRONGTYPE;
1856  }
1857 
1858  case SCIP_NLPPAR_VERBLEVEL:
1859  {
1860  SCIPerrorMessage("verbosity level parameter is of type int.\n");
1861  return SCIP_PARAMETERWRONGTYPE;
1862  }
1863 
1864  case SCIP_NLPPAR_FEASTOL:
1865  {
1866  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
1867  return SCIP_PARAMETERWRONGTYPE;
1868  }
1869 
1870  case SCIP_NLPPAR_RELOBJTOL:
1871  {
1872  SCIPerrorMessage("objective tolerance parameter is of type real.\n");
1873  return SCIP_PARAMETERWRONGTYPE;
1874  }
1875 
1876  case SCIP_NLPPAR_LOBJLIM:
1877  {
1878  SCIPerrorMessage("objective limit parameter is of type real.\n");
1879  return SCIP_PARAMETERWRONGTYPE;
1880  }
1881 
1882  case SCIP_NLPPAR_INFINITY:
1883  {
1884  SCIPerrorMessage("infinity parameter is of type real.\n");
1885  return SCIP_PARAMETERWRONGTYPE;
1886  }
1887 
1888  case SCIP_NLPPAR_ITLIM:
1889  {
1890  SCIPerrorMessage("iteration limit parameter is of type int.\n");
1891  return SCIP_PARAMETERWRONGTYPE;
1892  }
1893 
1894  case SCIP_NLPPAR_TILIM:
1895  {
1896  SCIPerrorMessage("time limit parameter is of type real.\n");
1897  return SCIP_PARAMETERWRONGTYPE;
1898  }
1899 
1900  case SCIP_NLPPAR_OPTFILE:
1901  {
1902  if( !problem->optfile.empty() )
1903  *sval = problem->optfile.c_str();
1904  else
1905  *sval = NULL;
1906  return SCIP_OKAY;
1907  }
1908 
1909  case SCIP_NLPPAR_FASTFAIL:
1910  {
1911  SCIPerrorMessage("fastfail parameter is of type int.\n");
1912  return SCIP_PARAMETERWRONGTYPE;
1913  }
1914 
1915  default:
1916  {
1917  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
1918  return SCIP_PARAMETERUNKNOWN;
1919  }
1920  }
1921 
1922  return SCIP_OKAY;
1923 }
1924 
1925 /** sets string parameter of NLP
1926  *
1927  * input:
1928  * - nlpi NLP interface structure
1929  * - problem datastructure for problem instance
1930  * - type parameter number
1931  * - sval parameter value
1932  */
1933 static
1934 SCIP_DECL_NLPISETSTRINGPAR( nlpiSetStringParIpopt )
1935 {
1936  assert(nlpi != NULL);
1937  assert(problem != NULL);
1938 
1939  switch( type )
1940  {
1942  {
1943  SCIPerrorMessage("from scratch parameter is of type int.\n");
1944  return SCIP_PARAMETERWRONGTYPE;
1945  }
1946 
1947  case SCIP_NLPPAR_VERBLEVEL:
1948  {
1949  SCIPerrorMessage("verbosity level parameter is of type int.\n");
1950  return SCIP_PARAMETERWRONGTYPE;
1951  }
1952 
1953  case SCIP_NLPPAR_FEASTOL:
1954  {
1955  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
1956  return SCIP_PARAMETERWRONGTYPE;
1957  }
1958 
1959  case SCIP_NLPPAR_RELOBJTOL:
1960  {
1961  SCIPerrorMessage("objective tolerance parameter is of type real.\n");
1962  return SCIP_PARAMETERWRONGTYPE;
1963  }
1964 
1965  case SCIP_NLPPAR_LOBJLIM:
1966  {
1967  SCIPerrorMessage("objective limit parameter is of type real.\n");
1968  return SCIP_PARAMETERWRONGTYPE;
1969  }
1970 
1971  case SCIP_NLPPAR_INFINITY:
1972  {
1973  SCIPerrorMessage("infinity parameter is of type real.\n");
1974  return SCIP_PARAMETERWRONGTYPE;
1975  }
1976 
1977  case SCIP_NLPPAR_ITLIM:
1978  {
1979  SCIPerrorMessage("iteration limit parameter is of type int.\n");
1980  return SCIP_PARAMETERWRONGTYPE;
1981  }
1982 
1983  case SCIP_NLPPAR_TILIM:
1984  {
1985  SCIPerrorMessage("time limit parameter is of type real.\n");
1986  return SCIP_PARAMETERWRONGTYPE;
1987  }
1988 
1989  case SCIP_NLPPAR_OPTFILE:
1990  {
1991  if( sval != NULL )
1992  problem->optfile = sval;
1993  else
1994  problem->optfile.clear();
1995 
1996  if( problem->ipopt->Initialize(problem->optfile) != Solve_Succeeded )
1997  {
1998  SCIPerrorMessage("Error initializing Ipopt using optionfile \"%s\"\n", problem->optfile.c_str());
1999  return SCIP_ERROR;
2000  }
2001  problem->ipopt->Options()->GetBoolValue("store_intermediate", problem->storeintermediate, "");
2002  problem->firstrun = TRUE;
2003 
2004  return SCIP_OKAY;
2005  }
2006 
2007  case SCIP_NLPPAR_FASTFAIL:
2008  {
2009  SCIPerrorMessage("fastfail parameter is of type int.\n");
2010  return SCIP_PARAMETERWRONGTYPE;
2011  }
2012 
2013  default:
2014  {
2015  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
2016  return SCIP_PARAMETERUNKNOWN;
2017  }
2018  }
2019 
2020  return SCIP_OKAY;
2021 }
2022 
2023 /** sets message handler for message output
2024  *
2025  * input:
2026  * - nlpi NLP interface structure
2027  * - messagehdlr SCIP message handler, or NULL to suppress all output
2028  */
2029 static
2030 SCIP_DECL_NLPISETMESSAGEHDLR( nlpiSetMessageHdlrIpopt )
2031 {
2032  SCIP_NLPIDATA* nlpidata;
2033 
2034  assert(nlpi != NULL);
2035 
2036  nlpidata = SCIPnlpiGetData(nlpi);
2037  assert(nlpidata != NULL);
2038 
2039  nlpidata->messagehdlr = messagehdlr;
2040 
2041  return SCIP_OKAY; /*lint !e527*/
2042 } /*lint !e715*/
2043 
2044 /** create solver interface for Ipopt solver */
2046  BMS_BLKMEM* blkmem, /**< block memory data structure */
2047  SCIP_NLPI** nlpi /**< pointer to buffer for nlpi address */
2048  )
2049 {
2050  SCIP_NLPIDATA* nlpidata;
2051 
2052  assert(blkmem != NULL);
2053  assert(nlpi != NULL);
2054 
2055  SCIP_ALLOC( nlpidata = new SCIP_NLPIDATA(blkmem) );
2056 
2057  SCIP_CALL( SCIPnlpiCreate(nlpi,
2059  nlpiCopyIpopt, nlpiFreeIpopt, nlpiGetSolverPointerIpopt,
2060  nlpiCreateProblemIpopt, nlpiFreeProblemIpopt, nlpiGetProblemPointerIpopt,
2061  nlpiAddVarsIpopt, nlpiAddConstraintsIpopt, nlpiSetObjectiveIpopt,
2062  nlpiChgVarBoundsIpopt, nlpiChgConsSidesIpopt, nlpiDelVarSetIpopt, nlpiDelConstraintSetIpopt,
2063  nlpiChgLinearCoefsIpopt, nlpiChgQuadraticCoefsIpopt, nlpiChgExprtreeIpopt, nlpiChgNonlinCoefIpopt,
2064  nlpiChgObjConstantIpopt, nlpiSetInitialGuessIpopt, nlpiSolveIpopt, nlpiGetSolstatIpopt, nlpiGetTermstatIpopt,
2065  nlpiGetSolutionIpopt, nlpiGetStatisticsIpopt,
2066  nlpiGetWarmstartSizeIpopt, nlpiGetWarmstartMemoIpopt, nlpiSetWarmstartMemoIpopt,
2067  nlpiGetIntParIpopt, nlpiSetIntParIpopt, nlpiGetRealParIpopt, nlpiSetRealParIpopt,
2068  nlpiGetStringParIpopt, nlpiSetStringParIpopt, nlpiSetMessageHdlrIpopt,
2069  nlpidata) );
2070 
2071  return SCIP_OKAY;
2072 }
2073 
2074 /** gets string that identifies Ipopt (version number) */
2075 const char* SCIPgetSolverNameIpopt(void)
2076 {
2077  return "Ipopt " IPOPT_VERSION;
2078 }
2079 
2080 /** gets string that describes Ipopt */
2081 const char* SCIPgetSolverDescIpopt(void)
2082 {
2083  return "Interior Point Optimizer developed by A. Waechter et.al. (www.coin-or.org/Ipopt)";
2084 }
2085 
2086 /** returns whether Ipopt is available, i.e., whether it has been linked in */
2088 {
2089  return TRUE;
2090 }
2091 
2092 /** gives a pointer to the IpoptApplication object stored in Ipopt-NLPI's NLPI problem data structure */
2094  SCIP_NLPIPROBLEM* nlpiproblem /**< NLP problem of Ipopt-NLPI */
2095  )
2096 {
2097  assert(nlpiproblem != NULL);
2098 
2099  return (void*)GetRawPtr(nlpiproblem->ipopt);
2100 }
2101 
2102 /** gives a pointer to the NLPIORACLE object stored in Ipopt-NLPI's NLPI problem data structure */
2104  SCIP_NLPIPROBLEM* nlpiproblem /**< NLP problem of Ipopt-NLPI */
2105  )
2106 {
2107  assert(nlpiproblem != NULL);
2108 
2109  return nlpiproblem->oracle;
2110 }
2111 
2112 /** sets modified default settings that are used when setting up an Ipopt problem
2113  *
2114  * Do not forget to add a newline after the last option in optionsstring.
2115  */
2117  SCIP_NLPI* nlpi, /**< Ipopt NLP interface */
2118  const char* optionsstring, /**< string with options as in Ipopt options file */
2119  SCIP_Bool append /**< whether to append to modified default settings or to overwrite */
2120  )
2121 {
2122  SCIP_NLPIDATA* data;
2123 
2124  assert(nlpi != NULL);
2125 
2126  data = SCIPnlpiGetData(nlpi);
2127  assert(data != NULL);
2128 
2129  if( append )
2130  data->defoptions += optionsstring;
2131  else
2132  data->defoptions = optionsstring;
2133 }
2134 
2135 /** Method to return some info about the nlp */
2136 bool ScipNLP::get_nlp_info(
2137  Index& n, /**< place to store number of variables */
2138  Index& m, /**< place to store number of constraints */
2139  Index& nnz_jac_g, /**< place to store number of nonzeros in jacobian */
2140  Index& nnz_h_lag, /**< place to store number of nonzeros in hessian */
2141  IndexStyleEnum& index_style /**< place to store used index style (0-based or 1-based) */
2142  )
2143 {
2144  const int* offset;
2145  SCIP_RETCODE retcode;
2146 
2147  assert(nlpiproblem != NULL);
2148  assert(nlpiproblem->oracle != NULL);
2149 
2150  n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
2151  m = SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle);
2152 
2153  retcode = SCIPnlpiOracleGetJacobianSparsity(nlpiproblem->oracle, &offset, NULL);
2154  if( retcode != SCIP_OKAY )
2155  return false;
2156  assert(offset != NULL);
2157  nnz_jac_g = offset[m];
2158 
2159  if( !approxhessian )
2160  {
2161  retcode = SCIPnlpiOracleGetHessianLagSparsity(nlpiproblem->oracle, &offset, NULL);
2162  if( retcode != SCIP_OKAY )
2163  return false;
2164  assert(offset != NULL);
2165  nnz_h_lag = offset[n];
2166  }
2167  else
2168  {
2169  nnz_h_lag = 0;
2170  }
2171 
2172  index_style = TNLP::C_STYLE;
2173 
2174  return true;
2175 }
2176 
2177 /** Method to return the bounds for my problem */
2178 bool ScipNLP::get_bounds_info(
2179  Index n, /**< number of variables */
2180  Number* x_l, /**< buffer to store lower bounds on variables */
2181  Number* x_u, /**< buffer to store upper bounds on variables */
2182  Index m, /**< number of constraints */
2183  Number* g_l, /**< buffer to store lower bounds on constraints */
2184  Number* g_u /**< buffer to store lower bounds on constraints */
2185  )
2186 {
2187  assert(nlpiproblem != NULL);
2188  assert(nlpiproblem->oracle != NULL);
2189 
2190  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2191  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2192 
2193  assert(SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle) != NULL || n == 0);
2194  assert(SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle) != NULL || n == 0);
2195 
2196  BMScopyMemoryArray(x_l, SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle), n);
2197  BMScopyMemoryArray(x_u, SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle), n);
2198 #ifndef NDEBUG
2199  for( int i = 0; i < n; ++i )
2200  assert(x_l[i] <= x_u[i]);
2201 #endif
2202 
2203  for( int i = 0; i < m; ++i )
2204  {
2205  g_l[i] = SCIPnlpiOracleGetConstraintLhs(nlpiproblem->oracle, i);
2206  g_u[i] = SCIPnlpiOracleGetConstraintRhs(nlpiproblem->oracle, i);
2207  assert(g_l[i] <= g_u[i]);
2208  }
2209 
2210  return true;
2211 }
2212 
2213 /** Method to return the starting point for the algorithm */
2214 bool ScipNLP::get_starting_point(
2215  Index n, /**< number of variables */
2216  bool init_x, /**< whether initial values for primal values are requested */
2217  Number* x, /**< buffer to store initial primal values */
2218  bool init_z, /**< whether initial values for dual values of variable bounds are requested */
2219  Number* z_L, /**< buffer to store dual values for variable lower bounds */
2220  Number* z_U, /**< buffer to store dual values for variable upper bounds */
2221  Index m, /**< number of constraints */
2222  bool init_lambda, /**< whether initial values for dual values of constraints are required */
2223  Number* lambda /**< buffer to store dual values of constraints */
2224  )
2225 {
2226  assert(nlpiproblem != NULL);
2227  assert(nlpiproblem->oracle != NULL);
2228 
2229  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2230  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2231 
2232  if( init_x )
2233  {
2234  if( nlpiproblem->initguess )
2235  {
2236  BMScopyMemoryArray(x, nlpiproblem->initguess, n);
2237  }
2238  else
2239  {
2240  SCIP_Real lb, ub;
2241 
2242  SCIPdebugMessage("Ipopt started without initial primal values; make up starting guess by projecting 0 onto variable bounds\n");
2243 
2244  for( int i = 0; i < n; ++i )
2245  {
2246  lb = SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle)[i];
2247  ub = SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle)[i];
2248  if( lb > 0.0 )
2249  x[i] = SCIPrandomGetReal(randnumgen, lb, lb + MAXPERTURB*MIN(1.0, ub-lb));
2250  else if( ub < 0.0 )
2251  x[i] = SCIPrandomGetReal(randnumgen, ub - MAXPERTURB*MIN(1.0, ub-lb), ub);
2252  else
2253  x[i] = SCIPrandomGetReal(randnumgen, MAX(lb, -MAXPERTURB*MIN(1.0, ub-lb)), MIN(ub, MAXPERTURB*MIN(1.0, ub-lb)));
2254  }
2255  }
2256  }
2257  if( init_z || init_lambda )
2258  return false;
2259 
2260  return true;
2261 }
2262 
2263 /** Method to return the variables linearity. */
2264 bool ScipNLP::get_variables_linearity(
2265  Index n, /**< number of variables */
2266  LinearityType* var_types /**< buffer to store linearity types of variables */
2267  )
2268 {
2269  assert(nlpiproblem != NULL);
2270  assert(nlpiproblem->oracle != NULL);
2271 
2272  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2273 
2274  for( int i = 0; i < n; ++i )
2275  var_types[i] = (SCIPnlpiOracleGetVarDegree(nlpiproblem->oracle, i) <= 1 ? LINEAR : NON_LINEAR);
2276 
2277  return true;
2278 }
2279 
2280 /** Method to return the constraint linearity. */
2281 bool ScipNLP::get_constraints_linearity(
2282  Index m, /**< number of constraints */
2283  LinearityType* const_types /**< buffer to store linearity types of constraints */
2284  )
2285 {
2286  int i;
2287 
2288  assert(nlpiproblem != NULL);
2289  assert(nlpiproblem->oracle != NULL);
2290 
2291  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2292 
2293  for( i = 0; i < m; ++i )
2294  const_types[i] = (SCIPnlpiOracleGetConstraintDegree(nlpiproblem->oracle, i) <= 1 ? LINEAR : NON_LINEAR);
2295 
2296  return true;
2297 }
2298 
2299 /** Method to return the number of nonlinear variables. */
2300 Index ScipNLP::get_number_of_nonlinear_variables()
2301 {
2302  int count;
2303  int n;
2304 
2305  assert(nlpiproblem != NULL);
2306  assert(nlpiproblem->oracle != NULL);
2307 
2308  n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
2309 
2310  count = 0;
2311  for( int i = 0; i < n; ++i )
2312  if (SCIPnlpiOracleGetVarDegree(nlpiproblem->oracle, i) > 1)
2313  ++count;
2314 
2315  return count;
2316 }
2317 
2318 /** Method to return the indices of the nonlinear variables */
2319 bool ScipNLP::get_list_of_nonlinear_variables(
2320  Index num_nonlin_vars, /**< number of nonlinear variables */
2321  Index* pos_nonlin_vars /**< array to fill with indices of nonlinear variables */
2322  )
2323 {
2324  int count;
2325  int n;
2326 
2327  assert(nlpiproblem != NULL);
2328  assert(nlpiproblem->oracle != NULL);
2329 
2330  n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
2331 
2332  count = 0;
2333  for( int i = 0; i < n; ++i )
2334  if (SCIPnlpiOracleGetVarDegree(nlpiproblem->oracle, i) > 1)
2335  {
2336  assert(count < num_nonlin_vars);
2337  pos_nonlin_vars[count++] = i;
2338  }
2339 
2340  assert(count == num_nonlin_vars);
2341 
2342  return true;
2343 }
2344 
2345 /** Method to return metadata about variables and constraints */
2346 bool ScipNLP::get_var_con_metadata(
2347  Index n, /**< number of variables */
2348  StringMetaDataMapType& var_string_md, /**< variable meta data of string type */
2349  IntegerMetaDataMapType& var_integer_md,/**< variable meta data of integer type */
2350  NumericMetaDataMapType& var_numeric_md,/**< variable meta data of numeric type */
2351  Index m, /**< number of constraints */
2352  StringMetaDataMapType& con_string_md, /**< constraint meta data of string type */
2353  IntegerMetaDataMapType& con_integer_md,/**< constraint meta data of integer type */
2354  NumericMetaDataMapType& con_numeric_md /**< constraint meta data of numeric type */
2355  )
2356 {
2357  assert(nlpiproblem != NULL);
2358  assert(nlpiproblem->oracle != NULL);
2359  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2360  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2361 
2362  char** varnames = SCIPnlpiOracleGetVarNames(nlpiproblem->oracle);
2363  if( varnames != NULL )
2364  {
2365  std::vector<std::string>& varnamesvec(var_string_md["idx_names"]);
2366  varnamesvec.reserve(n);
2367  for( int i = 0; i < n; ++i )
2368  {
2369  if( varnames[i] != NULL )
2370  {
2371  varnamesvec.push_back(varnames[i]);
2372  }
2373  else
2374  {
2375  char buffer[20];
2376  sprintf(buffer, "nlpivar%8d", i);
2377  varnamesvec.push_back(buffer);
2378  }
2379  }
2380  }
2381 
2382  std::vector<std::string>& consnamesvec(con_string_md["idx_names"]);
2383  consnamesvec.reserve(m);
2384  for( int i = 0; i < m; ++i )
2385  {
2386  if( SCIPnlpiOracleGetConstraintName(nlpiproblem->oracle, i) != NULL )
2387  {
2388  consnamesvec.push_back(SCIPnlpiOracleGetConstraintName(nlpiproblem->oracle, i));
2389  }
2390  else
2391  {
2392  char buffer[20];
2393  sprintf(buffer, "nlpicons%8d", i);
2394  consnamesvec.push_back(buffer);
2395  }
2396  }
2397 
2398  return true;
2399 }
2400 
2401 /** Method to return the objective value */
2402 bool ScipNLP::eval_f(
2403  Index n, /**< number of variables */
2404  const Number* x, /**< point to evaluate */
2405  bool new_x, /**< whether some function evaluation method has been called for this point before */
2406  Number& obj_value /**< place to store objective function value */
2407  )
2408 {
2409  assert(nlpiproblem != NULL);
2410  assert(nlpiproblem->oracle != NULL);
2411 
2412  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2413 
2414  return SCIPnlpiOracleEvalObjectiveValue(nlpiproblem->oracle, x, &obj_value) == SCIP_OKAY;
2415 }
2416 
2417 /** Method to return the gradient of the objective */
2418 bool ScipNLP::eval_grad_f(
2419  Index n, /**< number of variables */
2420  const Number* x, /**< point to evaluate */
2421  bool new_x, /**< whether some function evaluation method has been called for this point before */
2422  Number* grad_f /**< buffer to store objective gradient */
2423  )
2424 {
2425  SCIP_Real dummy;
2426 
2427  assert(nlpiproblem != NULL);
2428  assert(nlpiproblem->oracle != NULL);
2429 
2430  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2431 
2432  return SCIPnlpiOracleEvalObjectiveGradient(nlpiproblem->oracle, x, TRUE, &dummy, grad_f) == SCIP_OKAY;
2433 }
2434 
2435 /** Method to return the constraint residuals */
2436 bool ScipNLP::eval_g(
2437  Index n, /**< number of variables */
2438  const Number* x, /**< point to evaluate */
2439  bool new_x, /**< whether some function evaluation method has been called for this point before */
2440  Index m, /**< number of constraints */
2441  Number* g /**< buffer to store constraint function values */
2442  )
2443 {
2444  assert(nlpiproblem != NULL);
2445  assert(nlpiproblem->oracle != NULL);
2446 
2447  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2448 
2449  return SCIPnlpiOracleEvalConstraintValues(nlpiproblem->oracle, x, g) == SCIP_OKAY;
2450 }
2451 
2452 /** Method to return:
2453  * 1) The structure of the jacobian (if "values" is NULL)
2454  * 2) The values of the jacobian (if "values" is not NULL)
2455  */
2456 bool ScipNLP::eval_jac_g(
2457  Index n, /**< number of variables */
2458  const Number* x, /**< point to evaluate */
2459  bool new_x, /**< whether some function evaluation method has been called for this point before */
2460  Index m, /**< number of constraints */
2461  Index nele_jac, /**< number of nonzero entries in jacobian */
2462  Index* iRow, /**< buffer to store row indices of nonzero jacobian entries, or NULL if values are requested */
2463  Index* jCol, /**< buffer to store column indices of nonzero jacobian entries, or NULL if values are requested */
2464  Number* values /**< buffer to store values of nonzero jacobian entries, or NULL if structure is requested */
2465  )
2466 {
2467  assert(nlpiproblem != NULL);
2468  assert(nlpiproblem->oracle != NULL);
2469 
2470  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2471  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2472 
2473  if( values == NULL )
2474  { /* Ipopt wants to know sparsity structure */
2475  const int* jacoffset;
2476  const int* jaccol;
2477  int j;
2478  int i;
2479 
2480  assert(iRow != NULL);
2481  assert(jCol != NULL);
2482 
2483  if( SCIPnlpiOracleGetJacobianSparsity(nlpiproblem->oracle, &jacoffset, &jaccol) != SCIP_OKAY )
2484  return false;
2485 
2486  assert(jacoffset[0] == 0);
2487  assert(jacoffset[m] == nele_jac);
2488  j = jacoffset[0];
2489  for( i = 0; i < m; ++i )
2490  for( ; j < jacoffset[i+1]; ++j )
2491  iRow[j] = i;
2492 
2493  BMScopyMemoryArray(jCol, jaccol, nele_jac);
2494  }
2495  else
2496  {
2497  if( SCIPnlpiOracleEvalJacobian(nlpiproblem->oracle, x, TRUE, NULL, values) != SCIP_OKAY )
2498  return false;
2499  }
2500 
2501  return true;
2502 }
2503 
2504 /** Method to return:
2505  * 1) The structure of the hessian of the lagrangian (if "values" is NULL)
2506  * 2) The values of the hessian of the lagrangian (if "values" is not NULL)
2507  */
2508 bool ScipNLP::eval_h(
2509  Index n, /**< number of variables */
2510  const Number* x, /**< point to evaluate */
2511  bool new_x, /**< whether some function evaluation method has been called for this point before */
2512  Number obj_factor, /**< weight for objective function */
2513  Index m, /**< number of constraints */
2514  const Number* lambda, /**< weights for constraint functions */
2515  bool new_lambda, /**< whether the hessian has been evaluated for these values of lambda before */
2516  Index nele_hess, /**< number of nonzero entries in hessian */
2517  Index* iRow, /**< buffer to store row indices of nonzero hessian entries, or NULL if values are requested */
2518  Index* jCol, /**< buffer to store column indices of nonzero hessian entries, or NULL if values are requested */
2519  Number* values /**< buffer to store values of nonzero hessian entries, or NULL if structure is requested */
2520  )
2521 {
2522  assert(nlpiproblem != NULL);
2523  assert(nlpiproblem->oracle != NULL);
2524 
2525  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2526  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2527 
2528  if( values == NULL )
2529  { /* Ipopt wants to know sparsity structure */
2530  const int* heslagoffset;
2531  const int* heslagcol;
2532  int j;
2533  int i;
2534 
2535  assert(iRow != NULL);
2536  assert(jCol != NULL);
2537 
2538  if( SCIPnlpiOracleGetHessianLagSparsity(nlpiproblem->oracle, &heslagoffset, &heslagcol) != SCIP_OKAY )
2539  return false;
2540 
2541  assert(heslagoffset[0] == 0);
2542  assert(heslagoffset[n] == nele_hess);
2543  j = heslagoffset[0];
2544  for( i = 0; i < n; ++i )
2545  for( ; j < heslagoffset[i+1]; ++j )
2546  iRow[j] = i;
2547 
2548  BMScopyMemoryArray(jCol, heslagcol, nele_hess);
2549  }
2550  else
2551  {
2552  if( SCIPnlpiOracleEvalHessianLag(nlpiproblem->oracle, x, TRUE, obj_factor, lambda, values) != SCIP_OKAY )
2553  return false;
2554  }
2555 
2556  return true;
2557 }
2558 
2559 /** Method called by the solver at each iteration.
2560  *
2561  * Checks whether Ctrl-C was hit.
2562  */
2563 bool ScipNLP::intermediate_callback(
2564  AlgorithmMode mode, /**< current mode of algorithm */
2565  Index iter, /**< current iteration number */
2566  Number obj_value, /**< current objective value */
2567  Number inf_pr, /**< current primal infeasibility */
2568  Number inf_du, /**< current dual infeasibility */
2569  Number mu, /**< current barrier parameter */
2570  Number d_norm, /**< current gradient norm */
2571  Number regularization_size,/**< current size of regularization */
2572  Number alpha_du, /**< current dual alpha */
2573  Number alpha_pr, /**< current primal alpha */
2574  Index ls_trials, /**< current number of linesearch trials */
2575  const IpoptData* ip_data, /**< pointer to Ipopt Data */
2576  IpoptCalculatedQuantities* ip_cq /**< pointer to current calculated quantities */
2577  )
2578 {
2579  if( nlpiproblem->storeintermediate && mode == RegularMode && inf_pr < nlpiproblem->lastsolinfeas )
2580  {
2581  Ipopt::TNLPAdapter* tnlp_adapter;
2582 
2583  tnlp_adapter = NULL;
2584  if( ip_cq != NULL )
2585  {
2586  Ipopt::OrigIpoptNLP* orignlp;
2587 
2588  orignlp = dynamic_cast<OrigIpoptNLP*>(GetRawPtr(ip_cq->GetIpoptNLP()));
2589  if( orignlp != NULL )
2590  tnlp_adapter = dynamic_cast<TNLPAdapter*>(GetRawPtr(orignlp->nlp()));
2591  }
2592 
2593  if( tnlp_adapter != NULL && ip_data != NULL && IsValid(ip_data->curr()) )
2594  {
2595  SCIPdebugMessage("update lastsol: inf_pr old = %g -> new = %g\n", nlpiproblem->lastsolinfeas, inf_pr);
2596 
2597  if( nlpiproblem->lastsolprimals == NULL )
2598  {
2599  assert(nlpiproblem->lastsoldualcons == NULL);
2600  assert(nlpiproblem->lastsoldualvarlb == NULL);
2601  assert(nlpiproblem->lastsoldualvarub == NULL);
2602  if( BMSallocMemoryArray(&nlpiproblem->lastsolprimals, SCIPnlpiOracleGetNVars(nlpiproblem->oracle)) == NULL ||
2603  BMSallocMemoryArray(&nlpiproblem->lastsoldualcons, SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle)) == NULL ||
2604  BMSallocMemoryArray(&nlpiproblem->lastsoldualvarlb, SCIPnlpiOracleGetNVars(nlpiproblem->oracle)) == NULL ||
2605  BMSallocMemoryArray(&nlpiproblem->lastsoldualvarub, SCIPnlpiOracleGetNVars(nlpiproblem->oracle)) == NULL )
2606  {
2607  SCIPerrorMessage("out-of-memory in ScipNLP::intermediate_callback()\n");
2608  return TRUE;
2609  }
2610  }
2611 
2612  assert(IsValid(ip_data->curr()->x()));
2613  tnlp_adapter->ResortX(*ip_data->curr()->x(), nlpiproblem->lastsolprimals);
2614  nlpiproblem->lastsolinfeas = inf_pr;
2615 
2616  assert(IsValid(ip_data->curr()->y_c()));
2617  assert(IsValid(ip_data->curr()->y_d()));
2618  tnlp_adapter->ResortG(*ip_data->curr()->y_c(), *ip_data->curr()->y_d(), nlpiproblem->lastsoldualcons);
2619 
2620  // need to clear arrays first because ResortBnds only sets values for non-fixed variables
2621  BMSclearMemoryArray(nlpiproblem->lastsoldualvarlb, SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2622  BMSclearMemoryArray(nlpiproblem->lastsoldualvarub, SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2623  assert(IsValid(ip_data->curr()->z_L()));
2624  assert(IsValid(ip_data->curr()->z_U()));
2625  tnlp_adapter->ResortBnds(*ip_data->curr()->z_L(), nlpiproblem->lastsoldualvarlb, *ip_data->curr()->z_U(), nlpiproblem->lastsoldualvarub);
2626 
2627  }
2628  }
2629 
2630  /* do convergence test if fastfail is enabled */
2631  if( nlpiproblem->fastfail )
2632  {
2633  int i;
2634 
2635  if( iter == 0 )
2636  {
2637  conv_lastrestoiter = -1;
2638  }
2639  else if( mode == RestorationPhaseMode )
2640  {
2641  conv_lastrestoiter = iter;
2642  }
2643  else if( conv_lastrestoiter == iter-1 )
2644  {
2645  /* just switched back from restoration mode, reset dual reduction targets */
2646  for( i = 0; i < convcheck_nchecks; ++i )
2647  conv_dutarget[i] = convcheck_minred[i] * inf_du;
2648  }
2649 
2650  if( iter == convcheck_startiter )
2651  {
2652  /* define initial targets and iteration limits */
2653  for( i = 0; i < convcheck_nchecks; ++i )
2654  {
2655  conv_prtarget[i] = convcheck_minred[i] * inf_pr;
2656  conv_dutarget[i] = convcheck_minred[i] * inf_du;
2657  conv_iterlim[i] = iter + convcheck_maxiter[i];
2658  }
2659  }
2660  else if( iter > convcheck_startiter )
2661  {
2662  /* check if we should stop */
2663  for( i = 0; i < convcheck_nchecks; ++i )
2664  {
2665  if( inf_pr <= conv_prtarget[i] )
2666  {
2667  /* sufficient reduction w.r.t. primal infeasibility target
2668  * reset target w.r.t. current infeasibilities
2669  */
2670  conv_prtarget[i] = convcheck_minred[i] * inf_pr;
2671  conv_dutarget[i] = convcheck_minred[i] * inf_du;
2672  conv_iterlim[i] = iter + convcheck_maxiter[i];
2673  }
2674  else if( iter >= conv_iterlim[i] )
2675  {
2676  /* we hit a limit, should we really stop? */
2677  SCIPdebugMessage("convcheck %d: inf_pr = %e > target %e; inf_du = %e target %e: ",
2678  i, inf_pr, conv_prtarget[i], inf_du, conv_dutarget[i]);
2679  if( mode == RegularMode && iter <= conv_lastrestoiter + convcheck_startiter )
2680  {
2681  /* if we returned from feasibility restoration recently, we allow some more iterations,
2682  * because Ipopt may go for optimality for some iterations, at the costs of infeasibility
2683  */
2684  SCIPdebugPrintf("continue, because restoration phase only %d iters ago\n", iter - conv_lastrestoiter);
2685  }
2686  else if( mode == RegularMode && inf_du <= conv_dutarget[i] && iter < conv_iterlim[i] + convcheck_maxiter[i] )
2687  {
2688  /* if dual reduction is sufficient, we allow for twice the number of iterations to reach primal infeas reduction */
2689  SCIPdebugPrintf("continue, because dual infeas. red. sufficient and only %d iters above limit\n", iter - conv_iterlim[i]);
2690  }
2691  else
2692  {
2693  SCIPdebugPrintf("abort\n");
2694  return false;
2695  }
2696  }
2697  }
2698  }
2699  }
2700 
2701  return (SCIPinterrupted() == FALSE);
2702 }
2703 
2704 /** This method is called when the algorithm is complete so the TNLP can store/write the solution. */
2705 void ScipNLP::finalize_solution(
2706  SolverReturn status, /**< solve and solution status */
2707  Index n, /**< number of variables */
2708  const Number* x, /**< primal solution values */
2709  const Number* z_L, /**< dual values of variable lower bounds */
2710  const Number* z_U, /**< dual values of variable upper bounds */
2711  Index m, /**< number of constraints */
2712  const Number* g, /**< values of constraints */
2713  const Number* lambda, /**< dual values of constraints */
2714  Number obj_value, /**< objective function value */
2715  const IpoptData* data, /**< pointer to Ipopt Data */
2716  IpoptCalculatedQuantities* cq /**< pointer to calculated quantities */
2717  )
2718 {
2719  assert(nlpiproblem != NULL);
2720  assert(nlpiproblem->oracle != NULL);
2721 
2722  assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2723  assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2724 
2725  bool check_feasibility = false; // whether we should check x for feasibility, if not NULL
2726  switch( status )
2727  {
2728  case SUCCESS:
2729  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_LOCOPT;
2730  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
2731  assert(x != NULL);
2732  break;
2733 
2734  case STOP_AT_ACCEPTABLE_POINT:
2735  /* if stop at acceptable point, then dual infeasibility can be arbitrary large, so claim only feasibility */
2736  case FEASIBLE_POINT_FOUND:
2737  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_FEASIBLE;
2738  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
2739  assert(x != NULL);
2740  break;
2741 
2742  case MAXITER_EXCEEDED:
2743  check_feasibility = true;
2744  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2745  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_ITLIM;
2746  break;
2747 
2748  case CPUTIME_EXCEEDED:
2749  check_feasibility = true;
2750  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2751  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_TILIM;
2752  break;
2753 
2754  case STOP_AT_TINY_STEP:
2755  case RESTORATION_FAILURE:
2756  case ERROR_IN_STEP_COMPUTATION:
2757  check_feasibility = true;
2758  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2759  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_NUMERR;
2760  break;
2761 
2762  case LOCAL_INFEASIBILITY:
2763  /* still check feasibility, since we let Ipopt solve with higher tolerance than actually required */
2764  check_feasibility = true;
2765  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_LOCINFEASIBLE;
2766  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
2767  break;
2768 
2769  case DIVERGING_ITERATES:
2770  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNBOUNDED;
2771  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
2772  break;
2773 
2774  case INVALID_NUMBER_DETECTED:
2775  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2776  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_EVALERR;
2777  break;
2778 
2779  case USER_REQUESTED_STOP:
2780  check_feasibility = true;
2781  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2782  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OKAY;
2783  break;
2784 
2785  case TOO_FEW_DEGREES_OF_FREEDOM:
2786  case INTERNAL_ERROR:
2787  case INVALID_OPTION:
2788  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2789  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OTHER;
2790  break;
2791 
2792  case OUT_OF_MEMORY:
2793  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2794  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_MEMERR;
2795  break;
2796 
2797  default:
2798  SCIPerrorMessage("Ipopt returned with unknown solution status %d\n", status);
2799  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2800  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_OTHER;
2801  break;
2802  }
2803 
2804  /* if Ipopt reports its solution as locally infeasible or we don't know feasibility, then report the intermediate point with lowest constraint violation, if available */
2805  if( (x == NULL || nlpiproblem->lastsolstat == SCIP_NLPSOLSTAT_LOCINFEASIBLE || nlpiproblem->lastsolstat == SCIP_NLPSOLSTAT_UNKNOWN) && nlpiproblem->lastsolinfeas != SCIP_INVALID )
2806  {
2807  /* if infeasibility of lastsol is not invalid, then lastsol values should exist */
2808  assert(nlpiproblem->lastsolprimals != NULL);
2809  assert(nlpiproblem->lastsoldualcons != NULL);
2810  assert(nlpiproblem->lastsoldualvarlb != NULL);
2811  assert(nlpiproblem->lastsoldualvarub != NULL);
2812 
2813  /* check if lastsol is feasible */
2814  Number constrvioltol;
2815  nlpiproblem->ipopt->Options()->GetNumericValue("acceptable_constr_viol_tol", constrvioltol, "");
2816  if( nlpiproblem->lastsolinfeas <= constrvioltol )
2817  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_FEASIBLE;
2818  else
2819  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2820 
2821  SCIPdebugMessage("drop Ipopt's final point and report intermediate locally %sfeasible solution with infeas %g instead (acceptable: %g)\n",
2822  nlpiproblem->lastsolstat == SCIP_NLPSOLSTAT_LOCINFEASIBLE ? "in" : "", nlpiproblem->lastsolinfeas, constrvioltol);
2823  }
2824  else
2825  {
2826  assert(x != NULL);
2827  assert(lambda != NULL);
2828  assert(z_L != NULL);
2829  assert(z_U != NULL);
2830 
2831  if( nlpiproblem->lastsolprimals == NULL )
2832  {
2833  assert(nlpiproblem->lastsoldualcons == NULL);
2834  assert(nlpiproblem->lastsoldualvarlb == NULL);
2835  assert(nlpiproblem->lastsoldualvarub == NULL);
2836  BMSallocMemoryArray(&nlpiproblem->lastsolprimals, n);
2837  BMSallocMemoryArray(&nlpiproblem->lastsoldualcons, m);
2838  BMSallocMemoryArray(&nlpiproblem->lastsoldualvarlb, n);
2839  BMSallocMemoryArray(&nlpiproblem->lastsoldualvarub, n);
2840 
2841  if( nlpiproblem->lastsolprimals == NULL || nlpiproblem->lastsoldualcons == NULL ||
2842  nlpiproblem->lastsoldualvarlb == NULL || nlpiproblem->lastsoldualvarub == NULL )
2843  {
2844  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2845  nlpiproblem->lasttermstat = SCIP_NLPTERMSTAT_MEMERR;
2846  return;
2847  }
2848  }
2849 
2850  BMScopyMemoryArray(nlpiproblem->lastsolprimals, x, n);
2851  BMScopyMemoryArray(nlpiproblem->lastsoldualcons, lambda, m);
2852  BMScopyMemoryArray(nlpiproblem->lastsoldualvarlb, z_L, n);
2853  BMScopyMemoryArray(nlpiproblem->lastsoldualvarub, z_U, n);
2854 
2855  if( check_feasibility && cq != NULL )
2856  {
2857  Number constrviol;
2858  Number constrvioltol;
2859 
2860  constrviol = cq->curr_constraint_violation();
2861 
2862  nlpiproblem->ipopt->Options()->GetNumericValue("acceptable_constr_viol_tol", constrvioltol, "");
2863  if( constrviol <= constrvioltol )
2864  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_FEASIBLE;
2865  else if( nlpiproblem->lastsolstat != SCIP_NLPSOLSTAT_LOCINFEASIBLE )
2866  nlpiproblem->lastsolstat = SCIP_NLPSOLSTAT_UNKNOWN;
2867  }
2868  }
2869 }
2870 
2871 /** Calls Lapacks Dsyev routine to compute eigenvalues and eigenvectors of a dense matrix.
2872  *
2873  * It's here, because we use Ipopt's interface to Lapack.
2874  */
2876  SCIP_Bool computeeigenvectors,/**< should also eigenvectors should be computed ? */
2877  int N, /**< dimension */
2878  SCIP_Real* a, /**< matrix data on input (size N*N); eigenvectors on output if computeeigenvectors == TRUE */
2879  SCIP_Real* w /**< buffer to store eigenvalues (size N) */
2880  )
2881 {
2882  int info;
2883 
2884  IpLapackDsyev(computeeigenvectors, N, a, N, w, info);
2885 
2886  if( info != 0 )
2887  {
2888  SCIPerrorMessage("There was an error when calling DSYEV. INFO = %d\n", info);
2889  return SCIP_ERROR;
2890  }
2891 
2892  return SCIP_OKAY;
2893 }
2894 
2895 /** solves a linear problem of the form Ax = b for a regular matrix 3*3 A */
2896 static
2898  SCIP_Real* A, /**< matrix data on input (size 3*3); filled column-wise */
2899  SCIP_Real* b, /**< right hand side vector (size 3) */
2900  SCIP_Real* x, /**< buffer to store solution (size 3) */
2901  SCIP_Bool* success /**< pointer to store if the solving routine was successful */
2902  )
2903 {
2904  SCIP_Real Acopy[9];
2905  SCIP_Real bcopy[3];
2906  int pivotcopy[3];
2907  const int N = 3;
2908  int info;
2909 
2910  assert(A != NULL);
2911  assert(b != NULL);
2912  assert(x != NULL);
2913  assert(success != NULL);
2914 
2915  BMScopyMemoryArray(Acopy, A, N*N);
2916  BMScopyMemoryArray(bcopy, b, N);
2917 
2918  /* compute the LU factorization */
2919  IpLapackDgetrf(N, Acopy, pivotcopy, N, info);
2920 
2921  if( info != 0 )
2922  {
2923  SCIPdebugMessage("There was an error when calling Dgetrf. INFO = %d\n", info);
2924  *success = FALSE;
2925  }
2926  else
2927  {
2928  *success = TRUE;
2929 
2930  /* solve linear problem */
2931  IpLapackDgetrs(N, 1, Acopy, N, pivotcopy, bcopy, N);
2932 
2933  /* copy the solution */
2934  BMScopyMemoryArray(x, bcopy, N);
2935  }
2936 
2937  return SCIP_OKAY;
2938 }
2939 
2940 /** solves a linear problem of the form Ax = b for a regular matrix A
2941  *
2942  * Calls Lapacks IpLapackDgetrf routine to calculate a LU factorization and uses this factorization to solve
2943  * the linear problem Ax = b.
2944  * It's here, because Ipopt is linked against Lapack.
2945  */
2947  int N, /**< dimension */
2948  SCIP_Real* A, /**< matrix data on input (size N*N); filled column-wise */
2949  SCIP_Real* b, /**< right hand side vector (size N) */
2950  SCIP_Real* x, /**< buffer to store solution (size N) */
2951  SCIP_Bool* success /**< pointer to store if the solving routine was successful */
2952  )
2953 {
2954  SCIP_Real* Acopy;
2955  SCIP_Real* bcopy;
2956  int* pivotcopy;
2957  int info;
2958 
2959  assert(N > 0);
2960  assert(A != NULL);
2961  assert(b != NULL);
2962  assert(x != NULL);
2963  assert(success != NULL);
2964 
2965  /* call SCIPsolveLinearProb3() for performance reasons */
2966  if( N == 3 )
2967  {
2968  SCIP_CALL( SCIPsolveLinearProb3(A, b, x, success) );
2969  return SCIP_OKAY;
2970  }
2971 
2972  Acopy = NULL;
2973  bcopy = NULL;
2974  pivotcopy = NULL;
2975 
2976  SCIP_ALLOC( BMSduplicateMemoryArray(&Acopy, A, N*N) );
2977  SCIP_ALLOC( BMSduplicateMemoryArray(&bcopy, b, N) );
2978  SCIP_ALLOC( BMSallocMemoryArray(&pivotcopy, N) );
2979 
2980  /* compute the LU factorization */
2981  IpLapackDgetrf(N, Acopy, pivotcopy, N, info);
2982 
2983  if( info != 0 )
2984  {
2985  SCIPdebugMessage("There was an error when calling Dgetrf. INFO = %d\n", info);
2986  *success = FALSE;
2987  }
2988  else
2989  {
2990  *success = TRUE;
2991 
2992  /* solve linear problem */
2993  IpLapackDgetrs(N, 1, Acopy, N, pivotcopy, bcopy, N);
2994 
2995  /* copy the solution */
2996  BMScopyMemoryArray(x, bcopy, N);
2997  }
2998 
2999  BMSfreeMemoryArray(&pivotcopy);
3000  BMSfreeMemoryArray(&bcopy);
3001  BMSfreeMemoryArray(&Acopy);
3002 
3003  return SCIP_OKAY;
3004 }
SCIP_RETCODE SCIPnlpiOracleEvalConstraintValues(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Real *convals)
Definition: nlpioracle.c:2438
enum SCIP_NlpTermStat SCIP_NLPTERMSTAT
Definition: type_nlpi.h:84
SCIP_Real * lastsoldualvarlb
Definition: nlpi_ipopt.cpp:153
const SCIP_Real * SCIPnlpiOracleGetVarUbs(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2228
#define BMSfreeMemoryArrayNull(ptr)
Definition: memory.h:140
methods to interpret (evaluate) an expression tree "fast"
SCIP_RETCODE SCIPnlpiOracleGetJacobianSparsity(SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
Definition: nlpioracle.c:2513
SCIP_EXPRINTCAPABILITY SCIPnlpiOracleGetEvalCapability(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2378
SmartPtr< ScipNLP > nlp
Definition: nlpi_ipopt.cpp:141
void * SCIPgetNlpiOracleIpopt(SCIP_NLPIPROBLEM *nlpiproblem)
SCIP_RETCODE SCIPnlpiOracleGetHessianLagSparsity(SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
Definition: nlpioracle.c:2793
SCIP_NLPSOLSTAT lastsolstat
Definition: nlpi_ipopt.cpp:149
static SCIP_DECL_NLPICREATEPROBLEM(nlpiCreateProblemIpopt)
Definition: nlpi_ipopt.cpp:575
int SCIPnlpiOracleGetNVars(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2198
#define infinity
Definition: gastrans.c:71
static SCIP_DECL_NLPIGETWARMSTARTMEMO(nlpiGetWarmstartMemoIpopt)
static SCIP_DECL_NLPIFREEPROBLEM(nlpiFreeProblemIpopt)
Definition: nlpi_ipopt.cpp:669
#define NLPI_NAME
Definition: nlpi_ipopt.cpp:69
struct SCIP_NlpiProblem SCIP_NLPIPROBLEM
Definition: type_nlpi.h:39
internal methods for NLPI solver interfaces
SCIP_RETCODE SCIPnlpiOracleSetObjective(SCIP_NLPIORACLE *oracle, const SCIP_Real constant, int nlin, const int *lininds, const SCIP_Real *linvals, int nquadelems, const SCIP_QUADELEM *quadelems, const int *exprvaridxs, const SCIP_EXPRTREE *exprtree)
Definition: nlpioracle.c:1610
SCIP_NLPIDATA * SCIPnlpiGetData(SCIP_NLPI *nlpi)
Definition: nlpi.c:735
methods to store an NLP and request function, gradient, and hessian values
void SCIPsetModifiedDefaultSettingsIpopt(SCIP_NLPI *nlpi, const char *optionsstring, SCIP_Bool append)
SCIP_Real lasttime
Definition: nlpi_ipopt.cpp:157
SCIP_NLPIORACLE * oracle
#define FALSE
Definition: def.h:73
SCIP_Real SCIPnlpiOracleGetConstraintLhs(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2279
static SCIP_DECL_NLPIFREE(nlpiFreeIpopt)
Definition: nlpi_ipopt.cpp:538
static SCIP_DECL_NLPIGETSOLUTION(nlpiGetSolutionIpopt)
SCIP_Bool SCIPisIpoptAvailableIpopt(void)
#define NLPI_DESC
Definition: nlpi_ipopt.cpp:70
SCIP_Real SCIPrandomGetReal(SCIP_RANDNUMGEN *randnumgen, SCIP_Real minrandval, SCIP_Real maxrandval)
Definition: misc.c:9967
#define FEASTOLFACTOR
Definition: nlpi_ipopt.cpp:81
#define TRUE
Definition: def.h:72
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:54
#define DEFAULT_PRINTLEVEL
Definition: nlpi_ipopt.cpp:76
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveValue(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Real *objval)
Definition: nlpioracle.c:2400
static SCIP_DECL_NLPIGETINTPAR(nlpiGetIntParIpopt)
#define BMSallocMemoryArray(ptr, num)
Definition: memory.h:115
static SCIP_DECL_NLPISETINITIALGUESS(nlpiSetInitialGuessIpopt)
void SCIPmessageVPrintInfo(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr, va_list ap)
Definition: message.c:599
#define SCIP_CALL_ABORT_QUIET(x)
Definition: def.h:338
static SCIP_DECL_NLPISETINTPAR(nlpiSetIntParIpopt)
#define SCIPdebugMessage
Definition: pub_message.h:87
SCIP_RETCODE SCIPnlpiSetRealPar(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_NLPPARAM type, SCIP_Real dval)
Definition: nlpi.c:672
unsigned int SCIP_EXPRINTCAPABILITY
static SCIP_DECL_NLPIGETPROBLEMPOINTER(nlpiGetProblemPointerIpopt)
Definition: nlpi_ipopt.cpp:701
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveGradient(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *objval, SCIP_Real *objgrad)
Definition: nlpioracle.c:2464
static SCIP_DECL_NLPISETOBJECTIVE(nlpiSetObjectiveIpopt)
Definition: nlpi_ipopt.cpp:805
SCIP_RETCODE SCIPnlpiOracleDelConsSet(SCIP_NLPIORACLE *oracle, int *delstats)
Definition: nlpioracle.c:1827
void * SCIPgetIpoptApplicationPointerIpopt(SCIP_NLPIPROBLEM *nlpiproblem)
static void setOpttol(SCIP_NLPIPROBLEM *nlpiproblem, SCIP_Real opttol)
Definition: nlpi_ipopt.cpp:488
SCIP_RETCODE SCIPnlpiOracleChgObjConstant(SCIP_NLPIORACLE *oracle, SCIP_Real objconstant)
Definition: nlpioracle.c:2182
static SCIP_DECL_NLPISETSTRINGPAR(nlpiSetStringParIpopt)
SCIP_VAR ** x
Definition: circlepacking.c:54
int SCIPnlpiOracleGetNConstraints(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2208
std::string optfile
Definition: nlpi_ipopt.cpp:142
SCIP_Real * lastsoldualvarub
Definition: nlpi_ipopt.cpp:154
SCIP_RETCODE SCIPnlpiSetMessageHdlr(SCIP_NLPI *nlpi, SCIP_MESSAGEHDLR *messagehdlr)
Definition: nlpi.c:722
#define DEFAULT_RANDSEED
Definition: nlpi_ipopt.cpp:83
SCIP_Real lastsolinfeas
Definition: nlpi_ipopt.cpp:155
#define DEFAULT_MAXITER
Definition: nlpi_ipopt.cpp:78
void SCIPmessageVPrintError(const char *formatstr, va_list ap)
Definition: message.c:795
SCIP_VAR * w
Definition: circlepacking.c:58
static SCIP_DECL_NLPIGETSTATISTICS(nlpiGetStatisticsIpopt)
int SCIPnlpiOracleGetConstraintDegree(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2320
#define SUCCESS
Definition: portab.h:41
#define BMSfreeMemoryArray(ptr)
Definition: memory.h:139
SCIP_RETCODE LapackDsyev(SCIP_Bool computeeigenvectors, int N, SCIP_Real *a, SCIP_Real *w)
SCIP_RETCODE SCIPnlpiOracleChgExprtree(SCIP_NLPIORACLE *oracle, int considx, const int *exprvaridxs, const SCIP_EXPRTREE *exprtree)
Definition: nlpioracle.c:2100
static SCIP_DECL_NLPICHGOBJCONSTANT(nlpiChgObjConstantIpopt)
SCIP_RETCODE SCIPnlpiOracleCreate(BMS_BLKMEM *blkmem, SCIP_NLPIORACLE **oracle)
Definition: nlpioracle.c:1328
SCIP_Bool firstrun
Definition: nlpi_ipopt.cpp:146
void SCIPmessagePrintError(const char *formatstr,...)
Definition: message.c:782
#define SCIPerrorMessage
Definition: pub_message.h:55
#define NLPI_PRIORITY
Definition: nlpi_ipopt.cpp:71
#define SCIPdebugPrintf
Definition: pub_message.h:90
const SCIP_Real * SCIPnlpiOracleGetVarLbs(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2218
SCIP_RETCODE SCIPnlpiOracleChgVarBounds(SCIP_NLPIORACLE *oracle, int nvars, const int *indices, const SCIP_Real *lbs, const SCIP_Real *ubs)
Definition: nlpioracle.c:1647
struct SCIP_NlpiData SCIP_NLPIDATA
Definition: type_nlpi.h:38
enum SCIP_NlpSolStat SCIP_NLPSOLSTAT
Definition: type_nlpi.h:69
static SCIP_DECL_NLPICHGQUADCOEFS(nlpiChgQuadraticCoefsIpopt)
Definition: nlpi_ipopt.cpp:977
static SCIP_DECL_NLPICHGNONLINCOEF(nlpiChgNonlinCoefIpopt)
SCIP_RETCODE SCIPnlpiOracleAddVars(SCIP_NLPIORACLE *oracle, int nvars, const SCIP_Real *lbs, const SCIP_Real *ubs, const char **varnames)
Definition: nlpioracle.c:1447
static SCIP_DECL_NLPIGETREALPAR(nlpiGetRealParIpopt)
const char * SCIPgetSolverDescIpopt(void)
void SCIPmessagePrintWarning(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr,...)
Definition: message.c:418
void SCIPnlpStatisticsSetNIterations(SCIP_NLPSTATISTICS *statistics, int niterations)
Definition: nlpi.c:838
internal miscellaneous methods
#define NULL
Definition: lpi_spx1.cpp:155
void SCIPrandomFree(SCIP_RANDNUMGEN **randnumgen, BMS_BLKMEM *blkmem)
Definition: misc.c:9929
static void invalidateSolution(SCIP_NLPIPROBLEM *problem)
Definition: nlpi_ipopt.cpp:433
#define SCIP_CALL(x)
Definition: def.h:364
SCIP_Real * lastsolprimals
Definition: nlpi_ipopt.cpp:151
void SCIPnlpStatisticsSetTotalTime(SCIP_NLPSTATISTICS *statistics, SCIP_Real totaltime)
Definition: nlpi.c:848
#define SCIP_DEFAULT_FEASTOL
Definition: def.h:171
static SCIP_DECL_NLPICHGVARBOUNDS(nlpiChgVarBoundsIpopt)
Definition: nlpi_ipopt.cpp:839
void SCIPmessagePrintInfo(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr,...)
Definition: message.c:585
static SCIP_DECL_NLPICOPY(nlpiCopyIpopt)
Definition: nlpi_ipopt.cpp:507
SCIP_Real SCIPnlpiOracleGetConstraintRhs(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2292
const char * SCIPgetSolverNameIpopt(void)
SCIP_Real SCIPnlpiOracleGetInfinity(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:1400
#define BMSduplicateMemoryArray(ptr, source, num)
Definition: memory.h:135
SCIP_RETCODE SCIPnlpiOracleFree(SCIP_NLPIORACLE **oracle)
Definition: nlpioracle.c:1355
static const int convcheck_nchecks
Definition: nlpi_ipopt.cpp:110
methods for catching the user CTRL-C interrupt
Ipopt NLP interface.
SmartPtr< IpoptApplication > ipopt
Definition: nlpi_ipopt.cpp:140
public data structures and miscellaneous methods
#define SCIP_EXPRINTCAPABILITY_HESSIAN
static SCIP_DECL_NLPISETWARMSTARTMEMO(nlpiSetWarmstartMemoIpopt)
#define SCIP_Bool
Definition: def.h:70
SCIP_Bool SCIPinterrupted(void)
Definition: interrupt.c:154
#define MAXPERTURB
Definition: nlpi_ipopt.cpp:80
SCIP_RETCODE SCIPnlpiOracleSetProblemName(SCIP_NLPIORACLE *oracle, const char *name)
Definition: nlpioracle.c:1412
static void setFeastol(SCIP_NLPIPROBLEM *nlpiproblem, SCIP_Real feastol)
Definition: nlpi_ipopt.cpp:457
SCIP_RETCODE SCIPnlpiOracleAddConstraints(SCIP_NLPIORACLE *oracle, int nconss, const SCIP_Real *lhss, const SCIP_Real *rhss, const int *nlininds, int *const *lininds, SCIP_Real *const *linvals, const int *nquadelems, SCIP_QUADELEM *const *quadelems, int *const *exprvaridxs, SCIP_EXPRTREE *const *exprtrees, const char **consnames)
Definition: nlpioracle.c:1532
static SCIP_DECL_NLPICHGLINEARCOEFS(nlpiChgLinearCoefsIpopt)
Definition: nlpi_ipopt.cpp:953
static const int convcheck_startiter
Definition: nlpi_ipopt.cpp:111
static SCIP_DECL_NLPIGETSOLSTAT(nlpiGetSolstatIpopt)
#define MAX(x, y)
Definition: tclique_def.h:83
SCIP_EXPORT SCIP_EXPRINTCAPABILITY SCIPexprintGetCapability(void)
static SCIP_DECL_NLPIGETTERMSTAT(nlpiGetTermstatIpopt)
SCIP_RETCODE SCIPnlpiOracleChgLinearCoefs(SCIP_NLPIORACLE *oracle, int considx, int nentries, const int *varidxs, const SCIP_Real *newcoefs)
Definition: nlpioracle.c:1905
#define BMScopyMemoryArray(ptr, source, num)
Definition: memory.h:126
static SCIP_DECL_NLPIADDCONSTRAINTS(nlpiAddConstraintsIpopt)
Definition: nlpi_ipopt.cpp:765
SCIP_RETCODE SCIPnlpiOracleChgQuadCoefs(SCIP_NLPIORACLE *oracle, int considx, int nquadelems, const SCIP_QUADELEM *quadelems)
Definition: nlpioracle.c:2002
SCIP_RETCODE SCIPnlpiOracleChgConsSides(SCIP_NLPIORACLE *oracle, int nconss, const int *indices, const SCIP_Real *lhss, const SCIP_Real *rhss)
Definition: nlpioracle.c:1683
SCIP_RETCODE SCIPnlpiCreate(SCIP_NLPI **nlpi, const char *name, const char *description, int priority, SCIP_DECL_NLPICOPY((*nlpicopy)), SCIP_DECL_NLPIFREE((*nlpifree)), SCIP_DECL_NLPIGETSOLVERPOINTER((*nlpigetsolverpointer)), SCIP_DECL_NLPICREATEPROBLEM((*nlpicreateproblem)), SCIP_DECL_NLPIFREEPROBLEM((*nlpifreeproblem)), SCIP_DECL_NLPIGETPROBLEMPOINTER((*nlpigetproblempointer)), SCIP_DECL_NLPIADDVARS((*nlpiaddvars)), SCIP_DECL_NLPIADDCONSTRAINTS((*nlpiaddconstraints)), SCIP_DECL_NLPISETOBJECTIVE((*nlpisetobjective)), SCIP_DECL_NLPICHGVARBOUNDS((*nlpichgvarbounds)), SCIP_DECL_NLPICHGCONSSIDES((*nlpichgconssides)), SCIP_DECL_NLPIDELVARSET((*nlpidelvarset)), SCIP_DECL_NLPIDELCONSSET((*nlpidelconsset)), SCIP_DECL_NLPICHGLINEARCOEFS((*nlpichglinearcoefs)), SCIP_DECL_NLPICHGQUADCOEFS((*nlpichgquadcoefs)), SCIP_DECL_NLPICHGEXPRTREE((*nlpichgexprtree)), SCIP_DECL_NLPICHGNONLINCOEF((*nlpichgnonlincoef)), SCIP_DECL_NLPICHGOBJCONSTANT((*nlpichgobjconstant)), SCIP_DECL_NLPISETINITIALGUESS((*nlpisetinitialguess)), SCIP_DECL_NLPISOLVE((*nlpisolve)), SCIP_DECL_NLPIGETSOLSTAT((*nlpigetsolstat)), SCIP_DECL_NLPIGETTERMSTAT((*nlpigettermstat)), SCIP_DECL_NLPIGETSOLUTION((*nlpigetsolution)), SCIP_DECL_NLPIGETSTATISTICS((*nlpigetstatistics)), SCIP_DECL_NLPIGETWARMSTARTSIZE((*nlpigetwarmstartsize)), SCIP_DECL_NLPIGETWARMSTARTMEMO((*nlpigetwarmstartmemo)), SCIP_DECL_NLPISETWARMSTARTMEMO((*nlpisetwarmstartmemo)), SCIP_DECL_NLPIGETINTPAR((*nlpigetintpar)), SCIP_DECL_NLPISETINTPAR((*nlpisetintpar)), SCIP_DECL_NLPIGETREALPAR((*nlpigetrealpar)), SCIP_DECL_NLPISETREALPAR((*nlpisetrealpar)), SCIP_DECL_NLPIGETSTRINGPAR((*nlpigetstringpar)), SCIP_DECL_NLPISETSTRINGPAR((*nlpisetstringpar)), SCIP_DECL_NLPISETMESSAGEHDLR((*nlpisetmessagehdlr)), SCIP_NLPIDATA *nlpidata)
Definition: nlpi.c:41
SCIP_RETCODE SCIPnlpiOracleDelVarSet(SCIP_NLPIORACLE *oracle, int *delstats)
Definition: nlpioracle.c:1717
static SCIP_DECL_NLPIGETSOLVERPOINTER(nlpiGetSolverPointerIpopt)
Definition: nlpi_ipopt.cpp:560
SCIP_VAR ** b
Definition: circlepacking.c:56
char * SCIPnlpiOracleGetConstraintName(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2305
#define SCIP_DEFAULT_INFINITY
Definition: def.h:168
#define SCIP_EXPRINTCAPABILITY_FUNCVALUE
#define SCIP_DEFAULT_DUALFEASTOL
Definition: def.h:174
SCIP_RETCODE SCIPnlpiOracleChgExprParam(SCIP_NLPIORACLE *oracle, int considx, int paramidx, SCIP_Real paramval)
Definition: nlpioracle.c:2157
int SCIPnlpiOracleGetVarDegree(SCIP_NLPIORACLE *oracle, int varidx)
Definition: nlpioracle.c:2250
static SCIP_DECL_NLPIGETWARMSTARTSIZE(nlpiGetWarmstartSizeIpopt)
SCIP_RETCODE SCIPnlpiOracleEvalHessianLag(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real objfactor, const SCIP_Real *lambda, SCIP_Real *hessian)
Definition: nlpioracle.c:2892
SCIP_RETCODE SCIPsolveLinearProb(int N, SCIP_Real *A, SCIP_Real *b, SCIP_Real *x, SCIP_Bool *success)
public methods for message output
SCIP_VAR * a
Definition: circlepacking.c:57
static SCIP_DECL_NLPIDELVARSET(nlpiDelVarSetIpopt)
Definition: nlpi_ipopt.cpp:898
#define SCIP_Real
Definition: def.h:163
SCIP_NLPTERMSTAT lasttermstat
Definition: nlpi_ipopt.cpp:150
static SCIP_DECL_NLPISETMESSAGEHDLR(nlpiSetMessageHdlrIpopt)
SCIP_Real * lastsoldualcons
Definition: nlpi_ipopt.cpp:152
SCIP_RETCODE SCIPrandomCreate(SCIP_RANDNUMGEN **randnumgen, BMS_BLKMEM *blkmem, unsigned int initialseed)
Definition: misc.c:9913
#define SCIP_INVALID
Definition: def.h:183
static SCIP_DECL_NLPISOLVE(nlpiSolveIpopt)
SCIP_RETCODE SCIPnlpiOracleSetInfinity(SCIP_NLPIORACLE *oracle, SCIP_Real infinity)
Definition: nlpioracle.c:1384
static SCIP_DECL_NLPIGETSTRINGPAR(nlpiGetStringParIpopt)
static SCIP_DECL_NLPIDELCONSSET(nlpiDelConstraintSetIpopt)
Definition: nlpi_ipopt.cpp:925
static SCIP_DECL_NLPISETREALPAR(nlpiSetRealParIpopt)
SCIP_RETCODE SCIPnlpiOracleEvalJacobian(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *convals, SCIP_Real *jacobi)
Definition: nlpioracle.c:2654
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:122
static const SCIP_Real convcheck_minred[convcheck_nchecks]
Definition: nlpi_ipopt.cpp:113
struct BMS_BlkMem BMS_BLKMEM
Definition: memory.h:429
char ** SCIPnlpiOracleGetVarNames(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2238
static SCIP_DECL_NLPICHGEXPRTREE(nlpiChgExprtreeIpopt)
Definition: nlpi_ipopt.cpp:998
#define SCIP_ALLOC(x)
Definition: def.h:375
#define SCIPABORT()
Definition: def.h:336
static SCIP_DECL_NLPICHGCONSSIDES(nlpiChgConsSidesIpopt)
Definition: nlpi_ipopt.cpp:873
static SCIP_RETCODE SCIPsolveLinearProb3(SCIP_Real *A, SCIP_Real *b, SCIP_Real *x, SCIP_Bool *success)
SCIP_RETCODE SCIPcreateNlpSolverIpopt(BMS_BLKMEM *blkmem, SCIP_NLPI **nlpi)
static SCIP_DECL_NLPIADDVARS(nlpiAddVarsIpopt)
Definition: nlpi_ipopt.cpp:720
#define SCIP_EXPRINTCAPABILITY_GRADIENT
int SCIPnlpiOracleGetMaxDegree(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2352
static const int convcheck_maxiter[convcheck_nchecks]
Definition: nlpi_ipopt.cpp:112