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