SCIP Doxygen Documentation
 
Loading...
Searching...
No Matches
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-2023 Zuse Institute Berlin (ZIB) */
7/* */
8/* Licensed under the Apache License, Version 2.0 (the "License"); */
9/* you may not use this file except in compliance with the License. */
10/* You may obtain a copy of the License at */
11/* */
12/* http://www.apache.org/licenses/LICENSE-2.0 */
13/* */
14/* Unless required by applicable law or agreed to in writing, software */
15/* distributed under the License is distributed on an "AS IS" BASIS, */
16/* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. */
17/* See the License for the specific language governing permissions and */
18/* limitations under the License. */
19/* */
20/* You should have received a copy of the Apache-2.0 license */
21/* along with SCIP; see the file LICENSE. If not visit scipopt.org. */
22/* */
23/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
24
25/**@file nlpi_ipopt.cpp
26 * @ingroup DEFPLUGINS_NLPI
27 * @brief Ipopt NLP interface
28 * @author Stefan Vigerske
29 * @author Benjamin Müller
30 *
31 * @todo if too few degrees of freedom, solve a slack-minimization problem instead?
32 * @todo automatically switch to Hessian approximation if Hessian is dense or slow? (only do so if opttol/solvertol is large?)
33 *
34 * This file can only be compiled if Ipopt is available.
35 * Otherwise, to resolve public functions, use nlpi_ipopt_dummy.c.
36 * Since the dummy code is C instead of C++, it has been moved into a separate file.
37 */
38
39/*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
40
41#include "scip/nlpi_ipopt.h"
42
43#include "scip/nlpioracle.h"
44#include "scip/exprinterpret.h"
45#include "scip/scip_nlpi.h"
47#include "scip/scip_mem.h"
48#include "scip/scip_message.h"
49#include "scip/scip_general.h"
50#include "scip/scip_numerics.h"
51#include "scip/scip_param.h"
52#include "scip/scip_solve.h"
53#include "scip/scip_copy.h"
54#include "scip/pub_misc.h"
55#include "scip/pub_paramset.h"
56
57#include <new> /* for std::bad_alloc */
58#include <sstream>
59#include <cstring>
60
61/* turn off some lint warnings for file */
62/*lint --e{1540,750,3701}*/
63
64#include "IpoptConfig.h"
65
66#if defined(__GNUC__) && IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
67#pragma GCC diagnostic ignored "-Wshadow"
68#endif
69#include "IpIpoptApplication.hpp"
70#include "IpIpoptCalculatedQuantities.hpp"
71#include "IpSolveStatistics.hpp"
72#include "IpJournalist.hpp"
73#include "IpIpoptData.hpp"
74#include "IpTNLPAdapter.hpp"
75#include "IpOrigIpoptNLP.hpp"
76#include "IpLapack.hpp"
77#if defined(__GNUC__) && IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
78#pragma GCC diagnostic warning "-Wshadow"
79#endif
80
81#if IPOPT_VERSION_MAJOR < 3 || (IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 12)
82#error "The Ipopt interface requires at least 3.12.0"
83#endif
84
85/* MUMPS that can be used by Ipopt is not threadsafe
86 * If we want SCIP to be threadsafe (SCIP_THREADSAFE), have std::mutex (C++11 or higher), and use Ipopt before 3.14,
87 * then we protect the call to Ipopt by a mutex if MUMPS is used as linear solver.
88 * Thus, we allow only one Ipopt run at a time.
89 * Ipopt 3.14 has this build-in to its MUMPS interface, so we won't have to take care of this.
90 */
91#if defined(SCIP_THREADSAFE) && __cplusplus >= 201103L && IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
92#define PROTECT_SOLVE_BY_MUTEX
93#include <mutex>
94static std::mutex solve_mutex; /*lint !e1756*/
95#endif
96
97using namespace Ipopt;
98
99#define NLPI_NAME "ipopt" /**< short concise name of solver */
100#define NLPI_DESC "Ipopt interface" /**< description of solver */
101#define NLPI_PRIORITY 1000 /**< priority */
102
103#define MAXPERTURB 0.01 /**< maximal perturbation of bounds in starting point heuristic */
104#define FEASTOLFACTOR 0.9 /**< factor for user-given feasibility tolerance to get feasibility tolerance that is actually passed to Ipopt */
105
106#define DEFAULT_RANDSEED 71 /**< initial random seed */
107
108// enable this to collect statistics on number of iterations and problem characteristics in csv-form in log
109// note that this overwrites given iterlimit
110// see https://git.zib.de/integer/scip/-/snippets/1213 for some script that evaluates the collected data
111// #define COLLECT_SOLVESTATS
112
113/* Convergence check (see ScipNLP::intermediate_callback)
114 *
115 * If the fastfail option is set to aggressive, then we stop Ipopt if the reduction in
116 * primal infeasibility is not sufficient for a consecutive number of iterations.
117 * With the parameters as given below, we require Ipopt to
118 * - not increase the primal infeasibility after 5 iterations
119 * - reduce the primal infeasibility by at least 50% within 10 iterations
120 * - reduce the primal infeasibility by at least 90% within 30 iterations
121 * The targets are updated once they are reached and the limit on allowed iterations to reach the new target is reset.
122 *
123 * In certain situations, it is allowed to exceed an iteration limit:
124 * - If we are in the first 10 (convcheck_startiter) iterations.
125 * - If we are within 10 (convcheck_startiter) iterations after the restoration phase ended.
126 * The reason for this is that during feasibility restoration phase Ipopt aims completely on
127 * reducing constraint violation, completely forgetting the objective function.
128 * When returning from feasibility restoration and considering the original objective again,
129 * it is unlikely that Ipopt will continue to decrease primal infeasibility, since it may now target on
130 * more on optimality again. Thus, we do not check convergence for a number of iterations.
131 * - If the target on dual infeasibility reduction has been achieved, we are below twice the iteration limit, and
132 * we are not in restoration mode.
133 * The reason for this is that if Ipopt makes good progress towards optimality,
134 * we want to allow some more iterations where primal infeasibility is not reduced.
135 * However, in restoration mode, dual infeasibility does not correspond to the original problem and
136 * the complete aim is to restore primal infeasibility.
137 */
138static const int convcheck_nchecks = 3; /**< number of convergence checks */
139static const int convcheck_startiter = 10; /**< iteration where to start convergence checking */
140static const int convcheck_maxiter[convcheck_nchecks] = { 5, 15, 30 }; /**< maximal number of iterations to achieve each convergence check */
141static const SCIP_Real convcheck_minred[convcheck_nchecks] = { 1.0, 0.5, 0.1 }; /**< minimal required infeasibility reduction in each convergence check */
142
143/// integer parameters of Ipopt to make available via SCIP parameters
144static const char* ipopt_int_params[] =
145 { "print_level" }; // print_level must be first
146
147/// string parameters of Ipopt to make available via SCIP parameters
148static const char* ipopt_string_params[] =
149 { "linear_solver",
150 "hsllib",
151 "pardisolib",
152 "linear_system_scaling",
153 "nlp_scaling_method",
154 "mu_strategy",
155 "hessian_approximation"
156 };
157
158class ScipNLP;
159
160struct SCIP_NlpiData
161{
162public:
163 char* optfile; /**< Ipopt options file to read */
164 int print_level; /**< print_level set via nlpi/ipopt/print_level option */
165 SCIP_Real warm_start_push; /**< value to use for Ipopt's warm_start_bound_push/frac options */
166
167 /** constructor */
168 explicit SCIP_NlpiData()
169 : optfile(NULL), print_level(-1), warm_start_push(1e-9)
170 { }
171};
172
173struct SCIP_NlpiProblem
174{
175public:
176 SCIP_NLPIORACLE* oracle; /**< Oracle-helper to store and evaluate NLP */
177 SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
178
179 SmartPtr<IpoptApplication> ipopt; /**< Ipopt application */
180 SmartPtr<ScipNLP> nlp; /**< NLP in Ipopt form */
181
182 bool firstrun; /**< whether the next NLP solve will be the first one */
183 bool samestructure;/**< whether the NLP solved next will still have the same (Ipopt-internal) structure (same number of variables, constraints, bounds, and nonzero pattern) */
184
185 SCIP_NLPSOLSTAT solstat; /**< status of current solution (if any) */
186 SCIP_NLPTERMSTAT termstat; /**< termination status of last solve (if any) */
187 bool solprimalvalid;/**< whether primal solution values are available (solprimals has meaningful values) */
188 bool solprimalgiven;/**< whether primal solution values were set by caller */
189 bool soldualvalid; /**< whether dual solution values are available (soldual* have meaningful values) */
190 bool soldualgiven; /**< whether dual solution values were set by caller */
191 SCIP_Real* solprimals; /**< primal solution values, if available */
192 SCIP_Real* soldualcons; /**< dual solution values of constraints, if available */
193 SCIP_Real* soldualvarlb; /**< dual solution values of variable lower bounds, if available */
194 SCIP_Real* soldualvarub; /**< dual solution values of variable upper bounds, if available */
195 SCIP_Real solobjval; /**< objective function value in solution from last run */
196 SCIP_Real solconsviol; /**< constraint violation of primal solution, if available */
197 SCIP_Real solboundviol; /**< variable bound violation of primal solution, if available */
198 int lastniter; /**< number of iterations in last run */
199 SCIP_Real lasttime; /**< time spend in last run */
200
201 /** constructor */
211};
212
213/** TNLP implementation for SCIPs NLP */
214class ScipNLP : public TNLP
215{
216private:
217 SCIP_NLPIPROBLEM* nlpiproblem; /**< NLPI problem data */
218 SCIP* scip; /**< SCIP data structure */
219 SCIP_NLPPARAM param; /**< NLP solve parameters */
220
221 SCIP_Real conv_prtarget[convcheck_nchecks]; /**< target primal infeasibility for each convergence check */
222 SCIP_Real conv_dutarget[convcheck_nchecks]; /**< target dual infeasibility for each convergence check */
223 int conv_iterlim[convcheck_nchecks]; /**< iteration number where target primal infeasibility should to be achieved */
224 int conv_lastrestoiter; /**< last iteration number in restoration mode, or -1 if none */
225
226 unsigned int current_x; /**< unique number that identifies current iterate (x): incremented when Ipopt calls with new_x=true */
227 unsigned int last_f_eval_x; /**< the number of the iterate for which the objective was last evaluated (eval_f) */
228 unsigned int last_g_eval_x; /**< the number of the iterate for which the constraints were last evaluated (eval_g) */
229
230public:
231 bool approxhessian; /**< do we tell Ipopt to approximate the hessian? (may also be false if user set to approx. hessian via option file) */
232
233 // cppcheck-suppress uninitMemberVar
234 /** constructor */
235 ScipNLP(
236 SCIP_NLPIPROBLEM* nlpiproblem_ = NULL,/**< NLPI problem data */
237 SCIP* scip_ = NULL /**< SCIP data structure */
238 )
239 : nlpiproblem(nlpiproblem_), scip(scip_),
240 conv_lastrestoiter(-1),
241 current_x(1), last_f_eval_x(0), last_g_eval_x(0),
242 approxhessian(false)
243 {
244 assert(scip != NULL);
245 }
246
247 /** destructor */
248 ~ScipNLP()
249 { /*lint --e{1540}*/
250 }
251
252 /** initialize for new solve */
253 void initializeSolve(
254 SCIP_NLPIPROBLEM* nlpiproblem_, /**< NLPI problem */
255 const SCIP_NLPPARAM& nlpparam /**< NLP solve parameters */
256 )
257 {
259 nlpiproblem = nlpiproblem_;
260 param = nlpparam;
261
262 // since we are about to start a new solve, use this opportunity to reset the counts on x
263 current_x = 1;
264 last_f_eval_x = 0;
265 last_g_eval_x = 0;
266 }
267
268 /** Method to return some info about the nlp */
269 bool get_nlp_info(
270 Index& n, /**< place to store number of variables */
271 Index& m, /**< place to store number of constraints */
272 Index& nnz_jac_g, /**< place to store number of nonzeros in jacobian */
273 Index& nnz_h_lag, /**< place to store number of nonzeros in hessian */
274 IndexStyleEnum& index_style /**< place to store used index style (0-based or 1-based) */
275 );
276
277 /** Method to return the bounds for my problem */
278 bool get_bounds_info(
279 Index n, /**< number of variables */
280 Number* x_l, /**< buffer to store lower bounds on variables */
281 Number* x_u, /**< buffer to store upper bounds on variables */
282 Index m, /**< number of constraints */
283 Number* g_l, /**< buffer to store lower bounds on constraints */
284 Number* g_u /**< buffer to store lower bounds on constraints */
285 );
286
287 /** Method to return the starting point for the algorithm */
288 bool get_starting_point(
289 Index n, /**< number of variables */
290 bool init_x, /**< whether initial values for primal values are requested */
291 Number* x, /**< buffer to store initial primal values */
292 bool init_z, /**< whether initial values for dual values of variable bounds are requested */
293 Number* z_L, /**< buffer to store dual values for variable lower bounds */
294 Number* z_U, /**< buffer to store dual values for variable upper bounds */
295 Index m, /**< number of constraints */
296 bool init_lambda, /**< whether initial values for dual values of constraints are required */
297 Number* lambda /**< buffer to store dual values of constraints */
298 );
299
300 /** Method to return the number of nonlinear variables. */
301 Index get_number_of_nonlinear_variables();
302
303 /** Method to return the indices of the nonlinear variables */
304 bool get_list_of_nonlinear_variables(
305 Index num_nonlin_vars, /**< number of nonlinear variables */
306 Index* pos_nonlin_vars /**< array to fill with indices of nonlinear variables */
307 );
308
309 /** Method to return metadata about variables and constraints */
310 bool get_var_con_metadata(
311 Index n, /**< number of variables */
312 StringMetaDataMapType& var_string_md, /**< variable meta data of string type */
313 IntegerMetaDataMapType& var_integer_md,/**< variable meta data of integer type */
314 NumericMetaDataMapType& var_numeric_md,/**< variable meta data of numeric type */
315 Index m, /**< number of constraints */
316 StringMetaDataMapType& con_string_md, /**< constraint meta data of string type */
317 IntegerMetaDataMapType& con_integer_md,/**< constraint meta data of integer type */
318 NumericMetaDataMapType& con_numeric_md /**< constraint meta data of numeric type */
319 );
320
321 /** Method to return the objective value */
322 bool eval_f(
323 Index n, /**< number of variables */
324 const Number* x, /**< point to evaluate */
325 bool new_x, /**< whether some function evaluation method has been called for this point before */
326 Number& obj_value /**< place to store objective function value */
327 );
328
329 /** Method to return the gradient of the objective */
330 bool eval_grad_f(
331 Index n, /**< number of variables */
332 const Number* x, /**< point to evaluate */
333 bool new_x, /**< whether some function evaluation method has been called for this point before */
334 Number* grad_f /**< buffer to store objective gradient */
335 );
336
337 /** Method to return the constraint residuals */
338 bool eval_g(
339 Index n, /**< number of variables */
340 const Number* x, /**< point to evaluate */
341 bool new_x, /**< whether some function evaluation method has been called for this point before */
342 Index m, /**< number of constraints */
343 Number* g /**< buffer to store constraint function values */
344 );
345
346 /** Method to return:
347 * 1) The structure of the jacobian (if "values" is NULL)
348 * 2) The values of the jacobian (if "values" is not NULL)
349 */
350 bool eval_jac_g(
351 Index n, /**< number of variables */
352 const Number* x, /**< point to evaluate */
353 bool new_x, /**< whether some function evaluation method has been called for this point before */
354 Index m, /**< number of constraints */
355 Index nele_jac, /**< number of nonzero entries in jacobian */
356 Index* iRow, /**< buffer to store row indices of nonzero jacobian entries, or NULL if values
357 * are requested */
358 Index* jCol, /**< buffer to store column indices of nonzero jacobian entries, or NULL if values
359 * are requested */
360 Number* values /**< buffer to store values of nonzero jacobian entries, or NULL if structure is
361 * requested */
362 );
363
364 /** Method to return:
365 * 1) The structure of the hessian of the lagrangian (if "values" is NULL)
366 * 2) The values of the hessian of the lagrangian (if "values" is not NULL)
367 */
368 bool eval_h(
369 Index n, /**< number of variables */
370 const Number* x, /**< point to evaluate */
371 bool new_x, /**< whether some function evaluation method has been called for this point before */
372 Number obj_factor, /**< weight for objective function */
373 Index m, /**< number of constraints */
374 const Number* lambda, /**< weights for constraint functions */
375 bool new_lambda, /**< whether the hessian has been evaluated for these values of lambda before */
376 Index nele_hess, /**< number of nonzero entries in hessian */
377 Index* iRow, /**< buffer to store row indices of nonzero hessian entries, or NULL if values
378 * are requested */
379 Index* jCol, /**< buffer to store column indices of nonzero hessian entries, or NULL if values
380 * are requested */
381 Number* values /**< buffer to store values of nonzero hessian entries, or NULL if structure is requested */
382 );
383
384 /** Method called by the solver at each iteration.
385 *
386 * Checks whether Ctrl-C was hit.
387 */
388 bool intermediate_callback(
389 AlgorithmMode mode, /**< current mode of algorithm */
390 Index iter, /**< current iteration number */
391 Number obj_value, /**< current objective value */
392 Number inf_pr, /**< current primal infeasibility */
393 Number inf_du, /**< current dual infeasibility */
394 Number mu, /**< current barrier parameter */
395 Number d_norm, /**< current gradient norm */
396 Number regularization_size,/**< current size of regularization */
397 Number alpha_du, /**< current dual alpha */
398 Number alpha_pr, /**< current primal alpha */
399 Index ls_trials, /**< current number of linesearch trials */
400 const IpoptData* ip_data, /**< pointer to Ipopt Data */
401 IpoptCalculatedQuantities* ip_cq /**< pointer to current calculated quantities */
402 );
403
404 /** This method is called when the algorithm is complete so the TNLP can store/write the solution. */
405 void finalize_solution(
406 SolverReturn status, /**< solve and solution status */
407 Index n, /**< number of variables */
408 const Number* x, /**< primal solution values */
409 const Number* z_L, /**< dual values of variable lower bounds */
410 const Number* z_U, /**< dual values of variable upper bounds */
411 Index m, /**< number of constraints */
412 const Number* g, /**< values of constraints */
413 const Number* lambda, /**< dual values of constraints */
414 Number obj_value, /**< objective function value */
415 const IpoptData* data, /**< pointer to Ipopt Data */
416 IpoptCalculatedQuantities* cq /**< pointer to calculated quantities */
417 );
418};
419
420/** A particular Ipopt::Journal implementation that uses the SCIP message routines for output. */
421class ScipJournal : public Ipopt::Journal {
422private:
423 SCIP* scip; /**< SCIP data structure */
424
425public:
426 ScipJournal(
427 const char* name, /**< name of journal */
428 Ipopt::EJournalLevel default_level, /**< default verbosity level */
429 SCIP* scip_ /**< SCIP data structure */
430 )
431 : Ipopt::Journal(name, default_level),
432 scip(scip_)
433 { }
434
435 ~ScipJournal() { }
436
437protected:
438 /*lint -e{715}*/
439 void PrintImpl(
440 Ipopt::EJournalCategory category, /**< category of message */
441 Ipopt::EJournalLevel level, /**< verbosity level of message */
442 const char* str /**< message to print */
443 )
444 { /*lint --e{715} */
445 if( level == J_ERROR )
446 {
447 SCIPerrorMessage("%s", str);
448 }
449 else
450 {
451 SCIPinfoMessage(scip, NULL, "%s", str);
452 }
453 }
454
455 /*lint -e{715}*/
456 void PrintfImpl(
457 Ipopt::EJournalCategory category, /**< category of message */
458 Ipopt::EJournalLevel level, /**< verbosity level of message */
459 const char* pformat, /**< message printing format */
460 va_list ap /**< arguments of message */
461 )
462 { /*lint --e{715} */
463 if( level == J_ERROR )
464 {
467 }
468 else
469 {
471 }
472 }
473
474 void FlushBufferImpl() { }
475};
476
477/** sets status codes to mark that last NLP solve is no longer valid (usually because the NLP changed) */
478static
480 SCIP_NLPIPROBLEM* problem /**< data structure of problem */
481 )
482{
485 problem->solobjval = SCIP_INVALID;
486 problem->solconsviol = SCIP_INVALID;
487 problem->solboundviol = SCIP_INVALID;
488 problem->lastniter = -1;
489 problem->lasttime = -1.0;
490}
491
492/** sets solution values to be invalid and calls invalidateSolved() */
493static
495 SCIP_NLPIPROBLEM* problem /**< data structure of problem */
496 )
497{
498 assert(problem != NULL);
499
500 problem->solprimalvalid = false;
501 problem->solprimalgiven = false;
502 problem->soldualvalid = false;
503 problem->soldualgiven = false;
504
505 invalidateSolved(problem);
506}
507
508/** makes sure a starting point (initial guess) is available */
509static
511 SCIP* scip, /**< SCIP data structure */
512 SCIP_NLPIPROBLEM* problem, /**< data structure of problem */
513 SCIP_Bool& warmstart /**< whether a warmstart has been requested */
514 )
515{
516 SCIP_Real lb, ub;
517 int n;
518
519 assert(problem != NULL);
520
521 // disable warmstart if no primal or dual solution values are available
522 if( warmstart && (!problem->solprimalvalid || !problem->soldualvalid ))
523 {
524 SCIPdebugMsg(scip, "Disable warmstart as no primal or dual solution available.\n");
525 warmstart = false;
526 }
527
528 // continue below with making up a random primal starting point if
529 // the user did not set a starting point and warmstart is disabled (so the last solution shouldn't be used)
530 // (if warmstart, then due to the checks above we must now have valid primal and dual solution values)
531 if( problem->solprimalgiven || warmstart )
532 {
533 // so we must have a primal solution to start from
534 // if warmstart, then we also need to have a dual solution to start from
535 // if warmstart and primal solution is given by user, then also dual solution should have been given by user
536 assert(problem->solprimalvalid);
537 assert(problem->solprimals != NULL);
538 assert(!warmstart || !problem->solprimalgiven || problem->soldualgiven);
539 assert(!warmstart || problem->soldualcons != NULL);
540 assert(!warmstart || problem->soldualvarlb != NULL);
541 assert(!warmstart || problem->soldualvarub != NULL);
542 SCIPdebugMsg(scip, "Starting solution for %sstart available from %s.\n",
543 warmstart ? "warm" : "cold",
544 problem->solprimalgiven ? "user" : "previous solve");
545 return SCIP_OKAY;
546 }
547
548 SCIPdebugMsg(scip, "Starting solution for coldstart not available. Making up something by projecting 0 onto variable bounds and adding a random perturbation.\n");
549
550 n = SCIPnlpiOracleGetNVars(problem->oracle);
551
552 if( problem->randnumgen == NULL )
553 {
555 }
556
557 if( problem->solprimals == NULL )
558 {
560 }
561
562 for( int i = 0; i < n; ++i )
563 {
564 lb = SCIPnlpiOracleGetVarLbs(problem->oracle)[i];
565 ub = SCIPnlpiOracleGetVarUbs(problem->oracle)[i];
566 if( lb > 0.0 )
567 problem->solprimals[i] = SCIPrandomGetReal(problem->randnumgen, lb, lb + MAXPERTURB*MIN(1.0, ub-lb));
568 else if( ub < 0.0 )
569 problem->solprimals[i] = SCIPrandomGetReal(problem->randnumgen, ub - MAXPERTURB*MIN(1.0, ub-lb), ub);
570 else
571 problem->solprimals[i] = SCIPrandomGetReal(problem->randnumgen, MAX(lb, -MAXPERTURB*MIN(1.0, ub-lb)), MIN(ub, MAXPERTURB*MIN(1.0, ub-lb)));
572 }
573 problem->solprimalvalid = true;
574
575 return SCIP_OKAY;
576}
577
578/// pass NLP solve parameters to Ipopt
579static
581 SCIP* scip, /**< SCIP data structure */
582 SCIP_NLPIDATA* nlpidata, /**< NLPI data */
583 SCIP_NLPIPROBLEM* nlpiproblem, /**< NLP */
584 const SCIP_NLPPARAM param /**< solve parameters */
585 )
586{
587 assert(scip != NULL);
588 assert(nlpiproblem != NULL);
589
590 if( nlpidata->print_level < 0 ) // if nlpi/ipopt/print_level param has not been set
591 {
592 switch( param.verblevel )
593 {
594 case 0:
595 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("print_level", J_ERROR);
596 break;
597 case 1:
598 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("print_level", J_SUMMARY);
599 break;
600 case 2:
601 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("print_level", J_ITERSUMMARY);
602 break;
603 case 3:
604 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("print_level", J_DETAILED);
605 break;
606 default:
607 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("print_level", MIN(J_ITERSUMMARY + (param.verblevel-1), J_ALL));
608 break;
609 }
610 }
611
612#ifdef SCIP_DISABLED_CODE
613 if( param.iterlimit < 0 )
614 {
615 if( nlpidata->autoiterlim > 0 )
616 {
617 param.iterlimit = nlpidata->autoiterlim;
618 }
619 else
620 {
621 int nvars = 0; // number of variables, including slacks
622 int nnlvars = 0; // number of nonlinear variables
623 int nvarbnd = 0; // number of finite lower and upper bounds, including slacks
624 int nlincons = 0; // number of linear constraints
625 int nnlcons = 0; // number of nonlinear constraints
626 int jacnnz = 0; // number of nonzeros in Jacobian
627
628 int n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
629 int m = SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle);
630 const SCIP_Real* lbs = SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle);
631 const SCIP_Real* ubs = SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle);
632
633 for( int i = 0; i < n; ++i )
634 {
635 // skip fixed vars
636 if( SCIPisEQ(scip, lbs[i], ubs[i]) )
637 continue;
638
639 ++nvars;
640
641 if( SCIPnlpiOracleIsVarNonlinear(scip, nlpiproblem->oracle, i) )
642 ++nnlvars;
643
644 // every variable lower bound contributes a ln(x-lb) term in the barrier problem
645 if( !SCIPisInfinity(scip, -lbs[i]) )
646 ++nvarbnd;
647
648 // every variable upper bound contributes a ln(ub-x) term in the barrier problem
649 if( !SCIPisInfinity(scip, ubs[i]) )
650 ++nvarbnd;
651 }
652
653 for( int i = 0; i < m; ++i )
654 {
655 if( SCIPnlpiOracleIsConstraintNonlinear(nlpiproblem->oracle, i) )
656 ++nnlcons;
657 else
658 ++nlincons;
659
660 SCIP_Real lhs = SCIPnlpiOracleGetConstraintLhs(nlpiproblem->oracle, i);
661 SCIP_Real rhs = SCIPnlpiOracleGetConstraintRhs(nlpiproblem->oracle, i);
662
663 // every inequality contributes a slack variable
664 if( !SCIPisEQ(scip, lhs, rhs) )
665 ++nvars;
666
667 // every constraint lhs contributes a ln(slack-lhs) term in the barrier problem
668 if( !SCIPisInfinity(scip, -lhs) )
669 ++nvarbnd;
670 // every constraint rhs contributes a ln(rhs-slack) term in the barrier problem
671 if( !SCIPisInfinity(scip, rhs) )
672 ++nvarbnd;
673 }
674
675 const int* offset;
677 jacnnz = offset[m];
678
679 /* fitting data from NLP runs gave the following coefficients (see also !2634):
680 * intercept +40.2756726751
681 * jacnnz -0.0021259769
682 * jacnnz_sqrt +2.0121042012
683 * nlincons -0.0374801925
684 * nlincons_sqrt +2.9562232443
685 * nnlcons -0.0133039200
686 * nnlcons_sqrt -0.0412118434
687 * nnlvars -0.0702890379
688 * nnlvars_sqrt +7.0920920430
689 * nvarbnd +0.0183592749
690 * nvarbnd_sqrt -4.7218258847
691 * nvars +0.0112944627
692 * nvars_sqrt -0.8365873360
693 */
694 param.iterlimit = SCIPfloor(scip, 40.2756726751
695 -0.0021259769 * jacnnz +2.0121042012 * sqrt(jacnnz)
696 -0.0374801925 * nlincons +2.9562232443 * sqrt(nlincons)
697 -0.0133039200 * nnlcons -0.0412118434 * sqrt(nnlcons)
698 -0.0702890379 * nnlvars +7.0920920430 * sqrt(nnlvars)
699 +0.0183592749 * nvarbnd -4.7218258847 * sqrt(nvarbnd)
700 +0.0112944627 * nvars -0.8365873360 * sqrt(nvars));
701 SCIPdebugMsg(scip, "Iteration limit guess: %d\n", param.iterlimit);
702 /* but with all the negative coefficients, let's also ensure some minimal number of iterations */
703 if( param.iterlimit < 50 )
704 param.iterlimit = 50;
705 }
706 if( nlpidata->print_level >= J_SUMMARY || param.verblevel > 0 )
707 SCIPinfoMessage(scip, NULL, "Chosen iteration limit to be %d\n", param.iterlimit);
708 }
709#endif
710 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("max_iter", param.iterlimit);
711
712 (void) nlpiproblem->ipopt->Options()->SetNumericValue("constr_viol_tol", FEASTOLFACTOR * param.feastol);
713 (void) nlpiproblem->ipopt->Options()->SetNumericValue("acceptable_constr_viol_tol", FEASTOLFACTOR * param.feastol);
714
715 /* set optimality tolerance parameters in Ipopt
716 *
717 * Sets dual_inf_tol and compl_inf_tol to opttol and tol to solvertol.
718 * We leave acceptable_dual_inf_tol and acceptable_compl_inf_tol untouched for now, which means that if Ipopt has convergence problems, then
719 * it can stop with a solution that is still feasible, but essentially without a proof of local optimality.
720 * Note, that in this case we report only feasibility and not optimality of the solution (see ScipNLP::finalize_solution).
721 */
722 (void) nlpiproblem->ipopt->Options()->SetNumericValue("dual_inf_tol", param.opttol);
723 (void) nlpiproblem->ipopt->Options()->SetNumericValue("compl_inf_tol", param.opttol);
724 if( param.solvertol > 0.0 )
725 (void) nlpiproblem->ipopt->Options()->SetNumericValue("tol", param.solvertol);
726 else
727#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR > 14 || (IPOPT_VERSION_MINOR == 14 && IPOPT_VERSION_RELEASE >= 2)
728 (void) nlpiproblem->ipopt->Options()->UnsetValue("tol");
729#else
730 (void) nlpiproblem->ipopt->Options()->SetNumericValue("tol", 1e-8); // 1e-8 is the default
731#endif
732
733 /* Ipopt doesn't like a setting of exactly 0 for the max_*_time, so increase as little as possible in that case */
734#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR >= 14
735 (void) nlpiproblem->ipopt->Options()->SetNumericValue("max_wall_time", MAX(param.timelimit, DBL_MIN));
736#else
737 (void) nlpiproblem->ipopt->Options()->SetNumericValue("max_cpu_time", MAX(param.timelimit, DBL_MIN));
738#endif
739
740 // disable acceptable-point heuristic iff fastfail is completely off
741 // by default (fastfail=conservative), it seems useful to have Ipopt stop when it obviously doesn't make progress (like one of the NLPs in the bendersqp ctest)
743 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("acceptable_iter", 0);
744 else
745#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR > 14 || (IPOPT_VERSION_MINOR == 14 && IPOPT_VERSION_RELEASE >= 2)
746 (void) nlpiproblem->ipopt->Options()->UnsetValue("acceptable_iter");
747#else
748 (void) nlpiproblem->ipopt->Options()->SetIntegerValue("acceptable_iter", 15); // 15 is the default
749#endif
750
751 (void) nlpiproblem->ipopt->Options()->SetStringValue("expect_infeasible_problem", param.expectinfeas ? "yes" : "no");
752
753 if( !nlpiproblem->ipopt->Options()->SetStringValue("warm_start_init_point", param.warmstart ? "yes" : "no") && !param.warmstart )
754 {
755 // if we cannot disable warmstarts in Ipopt, then we have a big problem
756 SCIPerrorMessage("Failed to set Ipopt warm_start_init_point option to no.");
757 return SCIP_ERROR;
758 }
759
760 return SCIP_OKAY;
761}
762
763#ifdef COLLECT_SOLVESTATS
764/// writes out some solve status, number of iterations, time, and problem properties
765static
767 SCIP* scip,
769 SCIP_NLPIPROBLEM* problem,
771 )
772{
773 int nvars = 0; // number of variables, including slacks
774 int nslacks = 0; // number of slack variables
775 int nnlvars = 0; // number of nonlinear variables
776 int nvarlb = 0; // number of finite lower bounds, including slacks
777 int nvarub = 0; // number of finite upper bounds, including slacks
778 int nlincons = 0; // number of linear constraints
779 int nnlcons = 0; // number of nonlinear constraints
780 int objnl = 0; // number of nonlinear objectives
781 int jacnnz = 0; // number of nonzeros in Jacobian
782 int hesnnz = 0; // number of nonzeros in Hessian of Lagrangian
783 int linsys11nz = 0; // number of nonzeros in linear system (11) solved by Ipopt
784 int linsys13nz = 0; // number of nonzeros in linear system (13) solved by Ipopt
785
786 int n = SCIPnlpiOracleGetNVars(problem->oracle);
787 int m = SCIPnlpiOracleGetNConstraints(problem->oracle);
788
789 for( int i = 0; i < n; ++i )
790 {
791 SCIP_Real lb, ub;
792
793 lb = SCIPnlpiOracleGetVarLbs(problem->oracle)[i];
794 ub = SCIPnlpiOracleGetVarUbs(problem->oracle)[i];
795
796 // skip fixed vars
797 if( SCIPisEQ(scip, lb, ub) )
798 continue;
799
800 ++nvars;
801
802 if( SCIPnlpiOracleIsVarNonlinear(scip, problem->oracle, i) )
803 ++nnlvars;
804
805 // every variable lower bound contributes a ln(x-lb) term in the barrier problem
806 if( !SCIPisInfinity(scip, -lb) )
807 ++nvarlb;
808
809 // every variable upper bound contributes a ln(ub-x) term in the barrier problem
810 if( !SCIPisInfinity(scip, ub) )
811 ++nvarub;
812 }
813
814 for( int i = 0; i < m; ++i )
815 {
816 SCIP_Real lhs, rhs;
817
819 ++nnlcons;
820 else
821 ++nlincons;
822
823 lhs = SCIPnlpiOracleGetConstraintLhs(problem->oracle, i);
824 rhs = SCIPnlpiOracleGetConstraintRhs(problem->oracle, i);
825
826 // every inequality contributes a slack variable
827 if( !SCIPisEQ(scip, lhs, rhs) )
828 {
829 ++nvars;
830 ++nslacks;
831 }
832
833 // every constraint lhs contributes a ln(slack-lhs) term in the barrier problem
834 if( !SCIPisInfinity(scip, -lhs) )
835 ++nvarlb;
836 // every constraint rhs contributes a ln(rhs-slack) term in the barrier problem
837 if( !SCIPisInfinity(scip, rhs) )
838 ++nvarub;
839 }
840
842
843 const int* offset;
844 const int* col;
846 jacnnz = offset[m];
847
849 hesnnz = offset[n];
850
851 // number of nonzeros of matrix in linear system of barrier problem ((11) in Ipopt paper):
852 // off-diagonal elements of Hessian of Lagrangian (W_k without diagonal)
853 // + number of variables (diagonal of W_k + Sigma_k)
854 // + 2*(elements in Jacobian + number of slacks) (A_k)
855 for( int i = 0; i < n; ++i )
856 for( int j = offset[i]; j < offset[i+1]; ++j )
857 if( col[j] != i )
858 linsys11nz += 2; // off-diagonal element of Lagrangian, once for each triangle of matrix
859 linsys11nz += nvars; // number of variables
860 linsys11nz += 2 * (jacnnz + nslacks); // because each slack var contributes one entry to the Jacobian
861
862 // number of nonzeros of matrix in perturbed linear system of barrier problem ((13) in Ipopt paper):
863 // linsys11nz + number of constraints
865
866 SCIP_Real linsys11density = (SCIP_Real)linsys11nz / SQR((SCIP_Real)(nvars+m));
867 SCIP_Real linsys13density = (SCIP_Real)linsys13nz / SQR((SCIP_Real)(nvars+m));
868
869 bool expectinfeas;
870 problem->ipopt->Options()->GetBoolValue("expect_infeasible_problem", expectinfeas, "");
871
872 static bool firstwrite = true;
873 if( firstwrite )
874 {
875 printf("IPOPTSTAT status,iter,time,nvars,nnlvars,nvarlb,nvarub,nlincons,nnlcons,objnl,jacnnz,hesnnz,linsys11nz,linsys13nz,linsys11density,linsys13density,expectinfeas\n");
876 firstwrite = false;
877 }
878
879 printf("IPOPTSTAT %d,%d,%g,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%d,%f,%f,%d\n",
880 status, stats->IterationCount(), stats->TotalWallclockTime(),
882}
883#endif
884
885/** copy method of NLP interface (called when SCIP copies plugins) */
886static
893
894/** destructor of NLP interface to free nlpi data */
895static
897{
898 assert(nlpidata != NULL);
899
900 delete *nlpidata;
901 *nlpidata = NULL;
902
903 return SCIP_OKAY;
904}
905
906/** gets pointer for NLP solver to do dirty stuff */
907static
909{
910 if( problem == NULL )
911 return NULL;
912
913 return (void*)GetRawPtr(problem->ipopt);
914}
915
916/** creates a problem instance */
917static
919{
920 SCIP_NLPIDATA* data;
921
922 assert(nlpi != NULL);
923 assert(problem != NULL);
924
925 data = SCIPnlpiGetData(nlpi);
926 assert(data != NULL);
927
928 SCIP_ALLOC( *problem = new SCIP_NLPIPROBLEM ); /*lint !e774*/
929
930 SCIP_CALL( SCIPnlpiOracleCreate(scip, &(*problem)->oracle) );
931 SCIP_CALL( SCIPnlpiOracleSetProblemName(scip, (*problem)->oracle, name) );
932
933 try
934 {
935 /* initialize IPOPT without default journal */
936 (*problem)->ipopt = new IpoptApplication(false);
937
938 /* plugin our journal to get output through SCIP message handler */
939 SmartPtr<Journal> jrnl = new ScipJournal("console", J_ITERSUMMARY, scip);
940 jrnl->SetPrintLevel(J_DBG, J_NONE);
941 if( !(*problem)->ipopt->Jnlst()->AddJournal(jrnl) )
942 {
943 SCIPerrorMessage("Failed to register ScipJournal for IPOPT output.");
944 }
945
946 /* initialize Ipopt/SCIP NLP interface */
947 (*problem)->nlp = new ScipNLP(*problem, scip);
948 }
949 catch( const std::bad_alloc& )
950 {
951 SCIPerrorMessage("Not enough memory to initialize Ipopt.\n");
952 return SCIP_NOMEMORY;
953 }
954
955 for( size_t i = 0; i < sizeof(ipopt_string_params) / sizeof(const char*); ++i )
956 {
957 SCIP_PARAM* param;
959 char* paramval;
960
961 strcpy(paramname, "nlpi/" NLPI_NAME "/");
963 param = SCIPgetParam(scip, paramname);
964
965 // skip parameters that we didn't add to SCIP because they didn't exist in this build of Ipopt
966 if( param == NULL )
967 continue;
968
969 // if value wasn't left at the default, then pass to Ipopt and forbid overwriting
971 assert(paramval != NULL);
972 if( *paramval != '\0' )
973 (void) (*problem)->ipopt->Options()->SetStringValue(ipopt_string_params[i], paramval, false);
974 }
975
976 for( size_t i = 0; i < sizeof(ipopt_int_params) / sizeof(const char*); ++i )
977 {
978 SCIP_PARAM* param;
980 int paramval;
981
982 strcpy(paramname, "nlpi/" NLPI_NAME "/");
984 param = SCIPgetParam(scip, paramname);
985
986 // skip parameters that we didn't add to SCIP because they didn't exist in this build of Ipopt
987 if( param == NULL )
988 continue;
989
990 // if value wasn't left at the default, then pass to Ipopt and forbid overwriting
991 paramval = SCIPparamGetInt(param);
992 if( paramval != SCIPparamGetIntDefault(param) )
993 (void) (*problem)->ipopt->Options()->SetIntegerValue(ipopt_int_params[i], paramval, false);
994 }
995
996#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
997 /* Turn off bound relaxation for older Ipopt, as solutions may be out of bounds by more than constr_viol_tol.
998 * For Ipopt 3.14, bounds are relaxed by at most constr_viol_tol, so can leave bound_relax_factor at its default.
999 */
1000 (void) (*problem)->ipopt->Options()->SetNumericValue("bound_relax_factor", 0.0);
1001#endif
1002
1003 /* modify Ipopt's default settings to what we believe is appropriate */
1004 /* (*problem)->ipopt->Options()->SetStringValue("print_timing_statistics", "yes"); */
1005#ifdef SCIP_DEBUG
1006 (void) (*problem)->ipopt->Options()->SetStringValue("print_user_options", "yes");
1007#endif
1008 (void) (*problem)->ipopt->Options()->SetStringValue("sb", "yes");
1009 (void) (*problem)->ipopt->Options()->SetStringValueIfUnset("mu_strategy", "adaptive");
1010 (void) (*problem)->ipopt->Options()->SetIntegerValue("max_iter", INT_MAX);
1011 (void) (*problem)->ipopt->Options()->SetNumericValue("nlp_lower_bound_inf", -SCIPinfinity(scip), false);
1012 (void) (*problem)->ipopt->Options()->SetNumericValue("nlp_upper_bound_inf", SCIPinfinity(scip), false);
1013 (void) (*problem)->ipopt->Options()->SetNumericValue("diverging_iterates_tol", SCIPinfinity(scip), false);
1014
1015 /* when warmstarting, then reduce how much Ipopt modified the starting point */
1016 (void) (*problem)->ipopt->Options()->SetNumericValue("warm_start_bound_push", data->warm_start_push);
1017 (void) (*problem)->ipopt->Options()->SetNumericValue("warm_start_bound_frac", data->warm_start_push);
1018 (void) (*problem)->ipopt->Options()->SetNumericValue("warm_start_slack_bound_push", data->warm_start_push);
1019 (void) (*problem)->ipopt->Options()->SetNumericValue("warm_start_slack_bound_frac", data->warm_start_push);
1020 (void) (*problem)->ipopt->Options()->SetNumericValue("warm_start_mult_bound_push", data->warm_start_push);
1021
1022 /* apply user's given options file */
1023 assert(data->optfile != NULL);
1024 if( (*problem)->ipopt->Initialize(data->optfile) != Solve_Succeeded )
1025 {
1026 SCIPerrorMessage("Error during initialization of Ipopt using optionfile \"%s\"\n", data->optfile);
1027 return SCIP_ERROR;
1028 }
1029
1030 return SCIP_OKAY;
1031}
1032
1033/** free a problem instance */
1034static
1036{
1037 int n;
1038 int m;
1039
1040 assert(nlpi != NULL);
1041 assert(problem != NULL);
1042 assert(*problem != NULL);
1043 assert((*problem)->oracle != NULL);
1044
1045 n = SCIPnlpiOracleGetNVars((*problem)->oracle);
1046 m = SCIPnlpiOracleGetNConstraints((*problem)->oracle);
1047
1048 SCIPfreeBlockMemoryArrayNull(scip, &(*problem)->solprimals, n);
1049 SCIPfreeBlockMemoryArrayNull(scip, &(*problem)->soldualcons, m);
1050 SCIPfreeBlockMemoryArrayNull(scip, &(*problem)->soldualvarlb, n);
1051 SCIPfreeBlockMemoryArrayNull(scip, &(*problem)->soldualvarub, n);
1052
1053 SCIP_CALL( SCIPnlpiOracleFree(scip, &(*problem)->oracle) );
1054
1055 if( (*problem)->randnumgen != NULL )
1056 {
1057 SCIPfreeRandom(scip, &(*problem)->randnumgen);
1058 }
1059
1060 delete *problem;
1061 *problem = NULL;
1062
1063 return SCIP_OKAY;
1064}
1065
1066/** gets pointer to solver-internal problem instance to do dirty stuff */
1067static
1069{
1070 assert(nlpi != NULL);
1071 assert(problem != NULL);
1072
1073 return GetRawPtr(problem->nlp);
1074}
1075
1076/** add variables */
1077static
1079{
1080 int oldnvars;
1081
1082 assert(nlpi != NULL);
1083 assert(problem != NULL);
1084 assert(problem->oracle != NULL);
1085
1087
1091 invalidateSolution(problem);
1092
1093 SCIP_CALL( SCIPnlpiOracleAddVars(scip, problem->oracle, nvars, lbs, ubs, varnames) );
1094
1095 problem->samestructure = false;
1096
1097 return SCIP_OKAY;
1098}
1099
1100/** add constraints */
1101static
1103{
1104 int oldncons;
1105
1106 assert(nlpi != NULL);
1107 assert(problem != NULL);
1108 assert(problem->oracle != NULL);
1109
1111
1113 problem->soldualvalid = false;
1114 problem->soldualgiven = false;
1115
1117
1118 problem->samestructure = false;
1119
1120 return SCIP_OKAY;
1121}
1122
1123/** sets or overwrites objective, a minimization problem is expected */
1124static
1126{
1127 assert(nlpi != NULL);
1128 assert(problem != NULL);
1129 assert(problem->oracle != NULL);
1130
1131 /* We pass the objective gradient in dense form to Ipopt, so if the sparsity of that gradient changes, we do not change the structure of the problem inside Ipopt.
1132 * However, if the sparsity of the Hessian matrix of the objective changes, then the sparsity pattern of the Hessian of the Lagrangian may change.
1133 * Thus, set samestructure=false if the objective was and/or becomes nonlinear, but leave samestructure untouched if it was and stays linear.
1134 */
1135 if( expr != NULL || SCIPnlpiOracleIsConstraintNonlinear(problem->oracle, -1) )
1136 problem->samestructure = false;
1137
1138 SCIP_CALL( SCIPnlpiOracleSetObjective(scip, problem->oracle, constant, nlins, lininds, linvals, expr) );
1139
1140 /* keep solution as valid, but reset solve status and objective value */
1141 invalidateSolved(problem);
1142
1143 return SCIP_OKAY;
1144}
1145
1146/** change variable bounds */
1147static
1149{
1150 assert(nlpi != NULL);
1151 assert(problem != NULL);
1152 assert(problem->oracle != NULL);
1153
1154 /* Check whether the structure of the Ipopt internal NLP changes, if problem->samestructure at the moment.
1155 * We need to check whether variables become fixed or unfixed and whether bounds are added or removed.
1156 *
1157 * Project primal solution onto new bounds if currently valid.
1158 */
1159 if( problem->samestructure || problem->solprimalvalid )
1160 {
1161 for( int i = 0; i < nvars; ++i )
1162 {
1163 SCIP_Real oldlb;
1164 SCIP_Real oldub;
1165 oldlb = SCIPnlpiOracleGetVarLbs(problem->oracle)[indices[i]];
1166 oldub = SCIPnlpiOracleGetVarUbs(problem->oracle)[indices[i]];
1167
1168 if( (oldlb == oldub) != (lbs[i] == ubs[i]) ) /*lint !e777*/
1169 problem->samestructure = false;
1170 else if( SCIPisInfinity(scip, -oldlb) != SCIPisInfinity(scip, -lbs[i]) )
1171 problem->samestructure = false;
1172 else if( SCIPisInfinity(scip, oldub) != SCIPisInfinity(scip, ubs[i]) )
1173 problem->samestructure = false;
1174
1175 if( problem->solprimalvalid )
1176 {
1177 assert(problem->solprimals != NULL);
1178 problem->solprimals[i] = MIN(MAX(problem->solprimals[indices[i]], lbs[i]), ubs[i]);
1179 }
1180 }
1181 }
1182
1183 SCIP_CALL( SCIPnlpiOracleChgVarBounds(scip, problem->oracle, nvars, indices, lbs, ubs) );
1184
1185 invalidateSolved(problem);
1186
1187 return SCIP_OKAY;
1188}
1189
1190/** change constraint bounds */
1191static
1193{
1194 assert(nlpi != NULL);
1195 assert(problem != NULL);
1196 assert(problem->oracle != NULL);
1197
1198 /* Check whether the structure of the Ipopt internal NLP changes, if problem->samestructure at the moment.
1199 * We need to check whether constraints change from equality to inequality and whether sides are added or removed.
1200 */
1201 for( int i = 0; i < nconss && problem->samestructure; ++i )
1202 {
1203 SCIP_Real oldlhs;
1204 SCIP_Real oldrhs;
1205 oldlhs = SCIPnlpiOracleGetConstraintLhs(problem->oracle, indices[i]);
1206 oldrhs = SCIPnlpiOracleGetConstraintRhs(problem->oracle, indices[i]);
1207
1208 if( (oldlhs == oldrhs) != (lhss[i] == rhss[i]) ) /*lint !e777*/
1209 problem->samestructure = false;
1210 else if( SCIPisInfinity(scip, -oldlhs) != SCIPisInfinity(scip, -lhss[i]) )
1211 problem->samestructure = false;
1213 problem->samestructure = false;
1214 }
1215
1216 SCIP_CALL( SCIPnlpiOracleChgConsSides(scip, problem->oracle, nconss, indices, lhss, rhss) );
1217
1218 invalidateSolved(problem);
1219
1220 return SCIP_OKAY;
1221}
1222
1223/** delete a set of variables */
1224static
1226{
1227 int nvars;
1228
1229 assert(nlpi != NULL);
1230 assert(problem != NULL);
1231 assert(problem->oracle != NULL);
1233
1235
1237
1238 if( problem->solprimalvalid || problem->soldualvalid )
1239 {
1240 // update existing solution, if valid
1241 assert(!problem->solprimalvalid || problem->solprimals != NULL);
1242 assert(!problem->soldualvalid || problem->soldualvarlb != NULL);
1243 assert(!problem->soldualvalid || problem->soldualvarub != NULL);
1244
1245 int i;
1246 for( i = 0; i < dstatssize; ++i )
1247 {
1248 if( dstats[i] != -1 )
1249 {
1250 assert(dstats[i] >= 0);
1251 assert(dstats[i] < nvars);
1252 if( problem->solprimals != NULL )
1253 problem->solprimals[dstats[i]] = problem->solprimals[i];
1254 if( problem->soldualvarlb != NULL )
1255 {
1256 assert(problem->soldualvarub != NULL);
1257 problem->soldualvarlb[dstats[i]] = problem->soldualvarlb[i];
1258 problem->soldualvarub[dstats[i]] = problem->soldualvarub[i];
1259 }
1260 }
1261 }
1262 }
1263
1264 /* resize solution point arrays */
1265 if( problem->solprimals != NULL )
1266 {
1268 }
1269 if( problem->soldualvarlb != NULL )
1270 {
1272 }
1273 if( problem->soldualvarub != NULL )
1274 {
1276 }
1277
1278 problem->samestructure = false;
1279
1280 invalidateSolved(problem);
1281
1282 return SCIP_OKAY;
1283}
1284
1285/** delete a set of constraints */
1286static
1288{
1289 int ncons;
1290
1291 assert(nlpi != NULL);
1292 assert(problem != NULL);
1293 assert(problem->oracle != NULL);
1295
1297
1299
1300 if( problem->soldualvalid )
1301 {
1302 // update existing dual solution
1303 assert(problem->soldualcons != NULL);
1304
1305 int i;
1306 for( i = 0; i < dstatssize; ++i )
1307 {
1308 if( dstats[i] != -1 )
1309 {
1310 assert(dstats[i] >= 0);
1311 assert(dstats[i] < ncons);
1312 problem->soldualcons[dstats[i]] = problem->soldualcons[i];
1313 }
1314 }
1315 }
1316
1317 /* resize dual solution point array */
1318 if( problem->soldualcons != NULL )
1319 {
1321 }
1322
1323 problem->samestructure = false;
1324
1325 invalidateSolved(problem);
1326
1327 return SCIP_OKAY;
1328}
1329
1330/** change one linear coefficient in a constraint or objective */
1331static
1333{
1334 assert(nlpi != NULL);
1335 assert(problem != NULL);
1336 assert(problem->oracle != NULL);
1337
1338 SCIP_CALL( SCIPnlpiOracleChgLinearCoefs(scip, problem->oracle, idx, nvals, varidxs, vals) );
1339
1340 problem->samestructure = false; // nonzero patterns may have changed; TODO SCIPnlpiOracleChgLinearCoefs() should let us know
1341 invalidateSolved(problem);
1342
1343 return SCIP_OKAY;
1344}
1345
1346/** replaces the expression tree of a constraint or objective */
1347static
1349{
1350 assert(nlpi != NULL);
1351 assert(problem != NULL);
1352 assert(problem->oracle != NULL);
1353
1354 SCIP_CALL( SCIPnlpiOracleChgExpr(scip, problem->oracle, idxcons, expr) );
1355
1356 problem->samestructure = false; // nonzero patterns may have changed
1357 invalidateSolved(problem);
1358
1359 return SCIP_OKAY;
1360}
1361
1362/** change the constant offset in the objective */
1363static
1365{
1366 SCIP_Real oldconstant;
1367
1368 assert(nlpi != NULL);
1369 assert(problem != NULL);
1370 assert(problem->oracle != NULL);
1371
1373
1375
1376 if( problem->solobjval != SCIP_INVALID )
1377 problem->solobjval += objconstant - oldconstant;
1378
1379 return SCIP_OKAY;
1380}
1381
1382/** sets initial guess for primal variables */
1383static
1385{
1386 int nvars;
1387
1388 assert(nlpi != NULL);
1389 assert(problem != NULL);
1390 assert(problem->oracle != NULL);
1391
1393
1394 if( primalvalues != NULL )
1395 {
1396 // copy primal solution
1397 SCIPdebugMsg(scip, "set initial guess primal values to user-given\n");
1398 if( problem->solprimals == NULL )
1399 {
1401 }
1402 BMScopyMemoryArray(problem->solprimals, primalvalues, nvars);
1403 problem->solprimalvalid = true;
1404 problem->solprimalgiven = true;
1405 }
1406 else
1407 {
1408 // invalid current primal solution (if any)
1409 if( problem->solprimalvalid )
1410 {
1411 SCIPdebugMsg(scip, "invalidate initial guess primal values on user-request\n");
1412 }
1413 problem->solprimalvalid = false;
1414 problem->solprimalgiven = false;
1415 }
1416
1417 if( consdualvalues != NULL && varlbdualvalues != NULL && varubdualvalues != NULL )
1418 {
1419 // copy dual solution, if completely given
1420 SCIPdebugMsg(scip, "set initial guess dual values to user-given\n");
1422 if( problem->soldualcons == NULL )
1423 {
1425 }
1426 BMScopyMemoryArray(problem->soldualcons, consdualvalues, ncons);
1427
1428 assert((problem->soldualvarlb == NULL) == (problem->soldualvarub == NULL));
1429 if( problem->soldualvarlb == NULL )
1430 {
1433 }
1434 BMScopyMemoryArray(problem->soldualvarlb, varlbdualvalues, nvars);
1435 BMScopyMemoryArray(problem->soldualvarub, varubdualvalues, nvars);
1436 problem->soldualvalid = true;
1437 problem->soldualgiven = true;
1438 }
1439 else
1440 {
1441 // invalid current dual solution (if any)
1442 if( problem->soldualvalid )
1443 {
1444 SCIPdebugMsg(scip, "invalidate initial guess dual values\n");
1445 }
1446 problem->soldualvalid = false;
1447 problem->soldualgiven = false;
1448 }
1449
1450 return SCIP_OKAY;
1451}
1452
1453/** tries to solve NLP */
1454static
1456{
1457 SCIP_NLPIDATA* nlpidata;
1459
1460 assert(nlpi != NULL);
1461 assert(problem != NULL);
1462 assert(problem->oracle != NULL);
1463
1464 assert(IsValid(problem->ipopt));
1465 assert(IsValid(problem->nlp));
1466
1467 nlpidata = SCIPnlpiGetData(nlpi);
1468 assert(nlpidata != NULL);
1469
1470 // print parameters if either nlpi/ipopt/print_level has been set high enough or solve called with verblevel>0
1471 if( nlpidata->print_level >= J_SUMMARY || param.verblevel > 0 )
1472 {
1473 SCIPinfoMessage(scip, NULL, "Ipopt solve for problem %s at subSCIP depth %d", SCIPnlpiOracleGetProblemName(problem->oracle), SCIPgetSubscipDepth(scip));
1474 SCIPinfoMessage(scip, NULL, " with parameters " SCIP_NLPPARAM_PRINT(param));
1475 }
1476
1478
1479 if( param.timelimit == 0.0 )
1480 {
1481 /* there is nothing we can do if we are not given any time */
1482 problem->lastniter = 0;
1483 problem->lasttime = 0.0;
1486
1487 return SCIP_OKAY;
1488 }
1489
1490#ifdef COLLECT_SOLVESTATS
1491 // if collecting statistics on how many iterations one may need for a certain Ipopt problem,
1492 // do not use a tight iteration limit (but use some limit to get finished)
1493 param.iterlimit = 1000;
1494#endif
1495
1496 // change status info to unsolved, just in case
1497 invalidateSolved(problem);
1498
1499 // ensure a starting point is available
1500 // also disables param.warmstart if no warmstart available
1501 SCIP_CALL( ensureStartingPoint(scip, problem, param.warmstart) );
1502
1503 // tell NLP that we are about to start a new solve
1504 problem->nlp->initializeSolve(problem, param);
1505
1506 // set Ipopt parameters
1507 SCIP_CALL( handleNlpParam(scip, nlpidata, problem, param) );
1508
1509 try
1510 {
1511#ifdef PROTECT_SOLVE_BY_MUTEX
1512 /* lock solve_mutex if Ipopt is going to use Mumps as linear solver
1513 * unlocking will happen in the destructor of guard, which is called when this block is left
1514 */
1515 std::unique_lock<std::mutex> guard(solve_mutex, std::defer_lock); /*lint !e{728}*/
1516 std::string linsolver;
1517 (void) problem->ipopt->Options()->GetStringValue("linear_solver", linsolver, "");
1518 if( linsolver == "mumps" )
1519 guard.lock();
1520#endif
1521
1522 if( problem->firstrun )
1523 {
1525
1527
1528 /* if the expression interpreter or some user expression do not support function values and gradients and Hessians,
1529 * change NLP parameters or give an error
1530 */
1532 {
1535 {
1536 SCIPerrorMessage("Do not have expression interpreter that can compute function values and gradients. Cannot solve NLP with Ipopt.\n");
1539 return SCIP_OKAY;
1540 }
1541
1542 /* enable Hessian approximation if we are nonquadratic and the expression interpreter or user expression do not support Hessians */
1543 if( !(cap & SCIP_EXPRINTCAPABILITY_HESSIAN) )
1544 {
1545 (void) problem->ipopt->Options()->SetStringValueIfUnset("hessian_approximation", "limited-memory");
1546 problem->nlp->approxhessian = true;
1547 }
1548 else
1549 problem->nlp->approxhessian = false;
1550 }
1551
1552#ifdef SCIP_DEBUG
1553 problem->ipopt->Options()->SetStringValue("derivative_test", problem->nlp->approxhessian ? "first-order" : "second-order");
1554#endif
1555
1556 status = problem->ipopt->OptimizeTNLP(GetRawPtr(problem->nlp));
1557 }
1558 else
1559 {
1560 // TODO to be strict, we should check whether the eval capability has been changed and the Hessian approximation needs to be enabled (in which case we should call OptimizeTNLP instead)
1561 problem->ipopt->Options()->SetStringValue("warm_start_same_structure", problem->samestructure ? "yes" : "no");
1562 status = problem->ipopt->ReOptimizeTNLP(GetRawPtr(problem->nlp));
1563 }
1564
1565 // catch the very bad status codes
1566 switch( status ) {
1567 // everything better than Not_Enough_Degrees_Of_Freedom is a non-serious error
1568 case Solve_Succeeded:
1572 case Diverging_Iterates:
1576 case Restoration_Failed:
1579#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR >= 14
1580 case Maximum_WallTime_Exceeded: // new in Ipopt 3.14
1581 // if Ipopt >= 3.14, finalize_solution should always have been called if we get these status codes
1582 // this should have left us with some solution (unless we ran out of memory in finalize_solution)
1585#endif
1586 problem->firstrun = false;
1587 problem->samestructure = true;
1588 break;
1589
1593 SCIPdebugMsg(scip, "NLP has too few degrees of freedom.\n");
1594 break;
1595
1597 SCIPdebugMsg(scip, "Ipopt failed because of an invalid number in function or derivative value\n");
1599 /* Ipopt may or may not have called finalize solution
1600 * if it didn't, then we should still have SCIP_NLPSOLSTAT_UNKNOWN as set in the invalidateSolved() call above
1601 * if it did, then finalize_solution will have set SCIP_NLPSOLSTAT_UNKNOWN or SCIP_NLPSOLSTAT_FEASIBLE
1602 */
1604 break;
1605
1609 SCIPerrorMessage("Ipopt returned with status \"Insufficient Memory\"\n");
1610 return SCIP_NOMEMORY;
1611
1612 // really bad ones that could be something very unexpected going wrong within Ipopt
1614 case Internal_Error:
1617 SCIPerrorMessage("Ipopt returned with application return status %d\n", status);
1618 break;
1619
1620 // the really bad ones that indicate rather a programming error
1622 case Invalid_Option:
1626 SCIPerrorMessage("Ipopt returned with application return status %d\n", status);
1627 return SCIP_ERROR;
1628 }
1629
1630#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
1631 SmartPtr<SolveStatistics> stats = problem->ipopt->Statistics();
1632 /* Ipopt does not provide access to the statistics if there was a serious error */
1633 if( IsValid(stats) )
1634 {
1635 problem->lastniter = stats->IterationCount();
1636 problem->lasttime = stats->TotalWallclockTime();
1637
1638#ifdef COLLECT_SOLVESTATS
1639 collectStatistic(scip, status, problem, stats);
1640#endif
1641 }
1642#else
1643 SmartPtr<IpoptData> ip_data = problem->ipopt->IpoptDataObject();
1644 /* I don't think that there is a situation where ip_data is NULL, but check here anyway */
1645 if( IsValid(ip_data) )
1646 {
1647 problem->lastniter = ip_data->iter_count();
1648 problem->lasttime = ip_data->TimingStats().OverallAlgorithm().TotalWallclockTime();
1649 }
1650#endif
1651 else
1652 {
1653 problem->lastniter = 0;
1654 problem->lasttime = 0.0;
1655 }
1656 }
1657 catch( IpoptException& except )
1658 {
1659 SCIPerrorMessage("Ipopt returned with exception: %s\n", except.Message().c_str());
1660 return SCIP_ERROR;
1661 }
1662
1663 return SCIP_OKAY;
1664}
1665
1666/** gives solution status */
1667static
1669{
1670 assert(nlpi != NULL);
1671 assert(problem != NULL);
1672
1673 return problem->solstat;
1674}
1675
1676/** gives termination reason */
1677static
1679{
1680 assert(nlpi != NULL);
1681 assert(problem != NULL);
1682
1683 return problem->termstat;
1684}
1685
1686/** gives primal and dual solution values */
1687static
1689{
1690 assert(nlpi != NULL);
1691 assert(problem != NULL);
1692
1693 if( primalvalues != NULL )
1694 *primalvalues = problem->solprimals;
1695
1696 if( consdualvalues != NULL )
1697 *consdualvalues = problem->soldualcons;
1698
1699 if( varlbdualvalues != NULL )
1700 *varlbdualvalues = problem->soldualvarlb;
1701
1702 if( varubdualvalues != NULL )
1703 *varubdualvalues = problem->soldualvarub;
1704
1705 if( objval != NULL )
1706 *objval = problem->solobjval;
1707
1708 return SCIP_OKAY;
1709}
1710
1711/** gives solve statistics */
1712static
1714{
1715 assert(nlpi != NULL);
1716 assert(problem != NULL);
1718
1719 statistics->niterations = problem->lastniter;
1720 statistics->totaltime = problem->lasttime;
1721 statistics->evaltime = SCIPnlpiOracleGetEvalTime(scip, problem->oracle);
1722 statistics->consviol = problem->solconsviol;
1723 statistics->boundviol = problem->solboundviol;
1724
1725 return SCIP_OKAY;
1726}
1727
1728/** create solver interface for Ipopt solver and includes it into SCIP, if Ipopt is available */
1730 SCIP* scip /**< SCIP data structure */
1731 )
1732{
1733 SCIP_NLPIDATA* nlpidata;
1734 SCIP_Bool advanced = FALSE;
1735
1736 assert(scip != NULL);
1737
1738 SCIP_ALLOC( nlpidata = new SCIP_NLPIDATA() ); /*lint !e774*/
1739
1749 nlpidata) );
1750
1752
1753 SCIP_CALL( SCIPaddStringParam(scip, "nlpi/" NLPI_NAME "/optfile", "name of Ipopt options file",
1754 &nlpidata->optfile, FALSE, "", NULL, NULL) );
1755
1756 SCIP_CALL( SCIPaddRealParam(scip, "nlpi/" NLPI_NAME "/warm_start_push", "amount (relative and absolute) by which starting point is moved away from bounds in warmstarts",
1757 &nlpidata->warm_start_push, FALSE, 1e-9, 0.0, 1.0, NULL, NULL) );
1758
1760 IpoptApplication::RegisterAllIpoptOptions(reg_options);
1761
1762 for( size_t i = 0; i < sizeof(ipopt_string_params) / sizeof(const char*); ++i )
1763 {
1765
1766 // skip options not available with this build of Ipopt
1767 if( !IsValid(option) )
1768 continue;
1769
1770 assert(option->Type() == OT_String);
1771
1772 // prefix parameter name with nlpi/ipopt
1773 std::string paramname("nlpi/" NLPI_NAME "/");
1774 paramname += option->Name();
1775
1776 // initialize description with short description from Ipopt
1777 std::stringstream descr;
1778 descr << option->ShortDescription();
1779
1780 // add valid values to description, if there are more than one
1781 // the only case where there are less than 2 valid strings should be when anything is valid (in which case there is one valid string with value "*")
1782 std::vector<RegisteredOption::string_entry> validvals = option->GetValidStrings();
1783 if( validvals.size() > 1 )
1784 {
1785 descr << " Valid values if not empty:";
1786 for( std::vector<RegisteredOption::string_entry>::iterator val = validvals.begin(); val != validvals.end(); ++val )
1787 descr << ' ' << val->value_;
1788 }
1789
1790#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR >= 14
1791 // since Ipopt 3.14, Ipopt options have an advanced flag
1792 advanced = option->Advanced();
1793#endif
1794
1795 // we use the empty string as default to recognize later whether the user set has set the option
1796 SCIP_CALL( SCIPaddStringParam(scip, paramname.c_str(), descr.str().c_str(), NULL, advanced, "", NULL, NULL) );
1797 }
1798
1799 for( size_t i = 0; i < sizeof(ipopt_int_params) / sizeof(const char*); ++i )
1800 {
1801 assert(i > 0 || strcmp(ipopt_int_params[0], "print_level") == 0); // we assume print_level at index 0
1802
1804
1805 // skip options not available with this build of Ipopt
1806 if( !IsValid(option) )
1807 continue;
1808
1809 assert(option->Type() == OT_Integer);
1810
1811 // prefix parameter name with nlpi/ipopt
1812 std::string paramname("nlpi/" NLPI_NAME "/");
1813 paramname += option->Name();
1814
1815 int lower = option->LowerInteger();
1816 int upper = option->UpperInteger();
1817
1818 // we use value lower-1 as signal that the option was not modified by the user
1819 // for that, we require a finite lower bound
1820 assert(lower > INT_MIN);
1821
1822 // initialize description with short description from Ipopt
1823 std::stringstream descr;
1824 descr << option->ShortDescription();
1825 descr << ' ' << (lower-1) << " to use NLPI or Ipopt default.";
1826
1827#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR >= 14
1828 // since Ipopt 3.14, Ipopt options have an advanced flag
1829 advanced = option->Advanced();
1830#endif
1831
1832 // we use the empty string as default to recognize later whether the user set has set the option
1833 SCIP_CALL( SCIPaddIntParam(scip, paramname.c_str(), descr.str().c_str(),
1834 i == 0 ? &nlpidata->print_level : NULL, advanced,
1835 lower-1, lower-1, upper, NULL, NULL) );
1836 }
1837
1838 return SCIP_OKAY;
1839} /*lint !e429 */
1840
1841/** gets string that identifies Ipopt (version number) */
1842const char* SCIPgetSolverNameIpopt(void)
1843{
1844 return "Ipopt " IPOPT_VERSION;
1845}
1846
1847/** gets string that describes Ipopt */
1848const char* SCIPgetSolverDescIpopt(void)
1849{
1850 return "Interior Point Optimizer developed by A. Waechter et.al. (github.com/coin-or/Ipopt)";
1851}
1852
1853/** returns whether Ipopt is available, i.e., whether it has been linked in */
1855{
1856 return TRUE;
1857}
1858
1859/** gives a pointer to the NLPIORACLE object stored in Ipopt-NLPI's NLPI problem data structure */
1861 SCIP_NLPIPROBLEM* nlpiproblem /**< NLP problem of Ipopt-NLPI */
1862 )
1863{
1864 assert(nlpiproblem != NULL);
1865
1866 return nlpiproblem->oracle;
1867}
1868
1869/** Method to return some info about the nlp */
1870bool ScipNLP::get_nlp_info(
1871 Index& n, /**< place to store number of variables */
1872 Index& m, /**< place to store number of constraints */
1873 Index& nnz_jac_g, /**< place to store number of nonzeros in jacobian */
1874 Index& nnz_h_lag, /**< place to store number of nonzeros in hessian */
1875 IndexStyleEnum& index_style /**< place to store used index style (0-based or 1-based) */
1876 )
1877{
1878 const int* offset;
1879 SCIP_RETCODE retcode;
1880
1881 assert(nlpiproblem != NULL);
1882 assert(nlpiproblem->oracle != NULL);
1883
1884 n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
1885 m = SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle);
1886
1887 retcode = SCIPnlpiOracleGetJacobianSparsity(scip, nlpiproblem->oracle, &offset, NULL);
1888 if( retcode != SCIP_OKAY )
1889 return false;
1890 assert(offset != NULL);
1891 nnz_jac_g = offset[m];
1892
1893 if( !approxhessian )
1894 {
1895 retcode = SCIPnlpiOracleGetHessianLagSparsity(scip, nlpiproblem->oracle, &offset, NULL);
1896 if( retcode != SCIP_OKAY )
1897 return false;
1898 assert(offset != NULL);
1899 nnz_h_lag = offset[n];
1900 }
1901 else
1902 {
1903 nnz_h_lag = 0;
1904 }
1905
1906 index_style = TNLP::C_STYLE;
1907
1908 return true;
1909}
1910
1911/** Method to return the bounds for my problem */
1912bool ScipNLP::get_bounds_info(
1913 Index n, /**< number of variables */
1914 Number* x_l, /**< buffer to store lower bounds on variables */
1915 Number* x_u, /**< buffer to store upper bounds on variables */
1916 Index m, /**< number of constraints */
1917 Number* g_l, /**< buffer to store lower bounds on constraints */
1918 Number* g_u /**< buffer to store lower bounds on constraints */
1919 )
1920{
1921 const int* varlincounts;
1922 const int* varnlcounts;
1923
1924 assert(nlpiproblem != NULL);
1925 assert(nlpiproblem->oracle != NULL);
1926
1927 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
1928 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
1929
1930 assert(SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle) != NULL || n == 0);
1931 assert(SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle) != NULL || n == 0);
1932
1933 BMScopyMemoryArray(x_l, SCIPnlpiOracleGetVarLbs(nlpiproblem->oracle), n);
1934 BMScopyMemoryArray(x_u, SCIPnlpiOracleGetVarUbs(nlpiproblem->oracle), n);
1935#ifndef NDEBUG
1936 for( int i = 0; i < n; ++i )
1937 assert(x_l[i] <= x_u[i]);
1938#endif
1939
1940 SCIPnlpiOracleGetVarCounts(scip, nlpiproblem->oracle, &varlincounts, &varnlcounts);
1941
1942 /* Ipopt performs better when unused variables do not appear, which we can achieve by fixing them,
1943 * since Ipopts TNLPAdapter will hide them from Ipopts NLP. In the dual solution, bound multipliers (z_L, z_U)
1944 * for these variables should have value 0.0 (they are set to -grad Lagrangian).
1945 */
1946 for( int i = 0; i < n; ++i )
1947 {
1948 if( varlincounts[i] == 0 && varnlcounts[i] == 0 )
1949 {
1950 SCIPdebugMsg(scip, "fix unused variable x%d [%g,%g] to 0.0 or bound\n", i, x_l[i], x_u[i]);
1951 assert(x_l[i] <= x_u[i]);
1952 x_l[i] = x_u[i] = MAX(MIN(x_u[i], 0.0), x_l[i]);
1953 }
1954 }
1955
1956 for( int i = 0; i < m; ++i )
1957 {
1958 g_l[i] = SCIPnlpiOracleGetConstraintLhs(nlpiproblem->oracle, i);
1959 g_u[i] = SCIPnlpiOracleGetConstraintRhs(nlpiproblem->oracle, i);
1960 assert(g_l[i] <= g_u[i]);
1961 }
1962
1963 return true;
1964}
1965
1966/** Method to return the starting point for the algorithm */ /*lint -e{715}*/
1967bool ScipNLP::get_starting_point(
1968 Index n, /**< number of variables */
1969 bool init_x, /**< whether initial values for primal values are requested */
1970 Number* x, /**< buffer to store initial primal values */
1971 bool init_z, /**< whether initial values for dual values of variable bounds are requested */
1972 Number* z_L, /**< buffer to store dual values for variable lower bounds */
1973 Number* z_U, /**< buffer to store dual values for variable upper bounds */
1974 Index m, /**< number of constraints */
1975 bool init_lambda, /**< whether initial values for dual values of constraints are required */
1976 Number* lambda /**< buffer to store dual values of constraints */
1977 )
1978{ /*lint --e{715} */
1979 assert(nlpiproblem != NULL);
1980 assert(nlpiproblem->oracle != NULL);
1981
1982 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
1983 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
1984
1985 if( init_x )
1986 {
1987 assert(nlpiproblem->solprimalvalid);
1988 assert(nlpiproblem->solprimals != NULL);
1989 BMScopyMemoryArray(x, nlpiproblem->solprimals, n);
1990 }
1991
1992 if( init_z )
1993 {
1994 assert(nlpiproblem->soldualvalid);
1995 assert(nlpiproblem->soldualvarlb != NULL);
1996 assert(nlpiproblem->soldualvarub != NULL);
1997 BMScopyMemoryArray(z_L, nlpiproblem->soldualvarlb, n);
1998 BMScopyMemoryArray(z_U, nlpiproblem->soldualvarub, n);
1999 }
2000
2001 if( init_lambda )
2002 {
2003 assert(nlpiproblem->soldualvalid);
2004 assert(nlpiproblem->soldualcons != NULL);
2005 BMScopyMemoryArray(lambda, nlpiproblem->soldualcons, m);
2006 }
2007
2008 return true;
2009}
2010
2011/** Method to return the number of nonlinear variables. */
2012Index ScipNLP::get_number_of_nonlinear_variables()
2013{
2014 int count;
2015 int n;
2016
2017 assert(nlpiproblem != NULL);
2018 assert(nlpiproblem->oracle != NULL);
2019
2020 n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
2021
2022 count = 0;
2023 for( int i = 0; i < n; ++i )
2024 if( SCIPnlpiOracleIsVarNonlinear(scip, nlpiproblem->oracle, i) )
2025 ++count;
2026
2027 return count;
2028}
2029
2030/** Method to return the indices of the nonlinear variables */
2031bool ScipNLP::get_list_of_nonlinear_variables(
2032 Index num_nonlin_vars, /**< number of nonlinear variables */
2033 Index* pos_nonlin_vars /**< array to fill with indices of nonlinear variables */
2034 )
2035{
2036 int count;
2037 int n;
2038
2039 assert(nlpiproblem != NULL);
2040 assert(nlpiproblem->oracle != NULL);
2041
2042 n = SCIPnlpiOracleGetNVars(nlpiproblem->oracle);
2043
2044 count = 0;
2045 for( int i = 0; i < n; ++i )
2046 {
2047 if( SCIPnlpiOracleIsVarNonlinear(scip, nlpiproblem->oracle, i) )
2048 {
2049 assert(count < num_nonlin_vars);
2050 pos_nonlin_vars[count++] = i;
2051 }
2052 }
2053
2054 assert(count == num_nonlin_vars);
2055
2056 return true;
2057}
2058
2059/** Method to return metadata about variables and constraints */ /*lint -e{715}*/
2060bool ScipNLP::get_var_con_metadata(
2061 Index n, /**< number of variables */
2062 StringMetaDataMapType& var_string_md, /**< variable meta data of string type */
2063 IntegerMetaDataMapType& var_integer_md,/**< variable meta data of integer type */
2064 NumericMetaDataMapType& var_numeric_md,/**< variable meta data of numeric type */
2065 Index m, /**< number of constraints */
2066 StringMetaDataMapType& con_string_md, /**< constraint meta data of string type */
2067 IntegerMetaDataMapType& con_integer_md,/**< constraint meta data of integer type */
2068 NumericMetaDataMapType& con_numeric_md /**< constraint meta data of numeric type */
2069 )
2070{ /*lint --e{715}*/
2071 assert(nlpiproblem != NULL);
2072 assert(nlpiproblem->oracle != NULL);
2073 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2074 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2075
2076 char** varnames = SCIPnlpiOracleGetVarNames(nlpiproblem->oracle);
2077 if( varnames != NULL )
2078 {
2079 std::vector<std::string>& varnamesvec(var_string_md["idx_names"]);
2080 varnamesvec.reserve((size_t)n);
2081 for( int i = 0; i < n; ++i )
2082 {
2083 if( varnames[i] != NULL )
2084 {
2085 varnamesvec.push_back(varnames[i]); /*lint !e3701*/
2086 }
2087 else
2088 {
2089 char buffer[20];
2090 (void) SCIPsnprintf(buffer, 20, "nlpivar%8d", i);
2091 varnamesvec.push_back(buffer);
2092 }
2093 }
2094 }
2095
2096 std::vector<std::string>& consnamesvec(con_string_md["idx_names"]);
2097 consnamesvec.reserve((size_t)m);
2098 for( int i = 0; i < m; ++i )
2099 {
2100 if( SCIPnlpiOracleGetConstraintName(nlpiproblem->oracle, i) != NULL )
2101 {
2102 consnamesvec.push_back(SCIPnlpiOracleGetConstraintName(nlpiproblem->oracle, i));
2103 }
2104 else
2105 {
2106 char buffer[20];
2107 (void) SCIPsnprintf(buffer, 20, "nlpicons%8d", i);
2108 consnamesvec.push_back(buffer); /*lint !e3701*/
2109 }
2110 }
2111
2112 return true;
2113}
2114
2115/** Method to return the objective value */ /*lint -e{715}*/
2116bool ScipNLP::eval_f(
2117 Index n, /**< number of variables */
2118 const Number* x, /**< point to evaluate */
2119 bool new_x, /**< whether some function evaluation method has been called for this point before */
2120 Number& obj_value /**< place to store objective function value */
2121 )
2122{ /*lint --e{715}*/
2123 assert(nlpiproblem != NULL);
2124 assert(nlpiproblem->oracle != NULL);
2125
2126 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2127
2128 if( new_x )
2129 ++current_x;
2130 last_f_eval_x = current_x;
2131
2132 return SCIPnlpiOracleEvalObjectiveValue(scip, nlpiproblem->oracle, x, &obj_value) == SCIP_OKAY;
2133}
2134
2135/** Method to return the gradient of the objective */ /*lint -e{715}*/
2136bool ScipNLP::eval_grad_f(
2137 Index n, /**< number of variables */
2138 const Number* x, /**< point to evaluate */
2139 bool new_x, /**< whether some function evaluation method has been called for this point before */
2140 Number* grad_f /**< buffer to store objective gradient */
2141 )
2142{ /*lint --e{715}*/
2143 SCIP_Real dummy;
2144
2145 assert(nlpiproblem != NULL);
2146 assert(nlpiproblem->oracle != NULL);
2147
2148 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2149
2150 if( new_x )
2151 ++current_x;
2152 else
2153 {
2154 // pass new_x = TRUE to objective gradient eval iff we have not evaluated the objective function at this point yet
2155 new_x = last_f_eval_x < current_x;
2156 }
2157 // if we evaluate the objective gradient with new_x = true, then this will also evaluate the objective function
2158 // (and if we do with new_x = false, then we already have last_f_eval_x == current_x anyway)
2159 last_f_eval_x = current_x;
2160
2161 return SCIPnlpiOracleEvalObjectiveGradient(scip, nlpiproblem->oracle, x, new_x, &dummy, grad_f) == SCIP_OKAY;
2162}
2163
2164/** Method to return the constraint residuals */ /*lint -e{715}*/
2165bool ScipNLP::eval_g(
2166 Index n, /**< number of variables */
2167 const Number* x, /**< point to evaluate */
2168 bool new_x, /**< whether some function evaluation method has been called for this point before */
2169 Index m, /**< number of constraints */
2170 Number* g /**< buffer to store constraint function values */
2171 )
2172{ /*lint --e{715}*/
2173 assert(nlpiproblem != NULL);
2174 assert(nlpiproblem->oracle != NULL);
2175
2176 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2177
2178 if( new_x )
2179 ++current_x;
2180 last_g_eval_x = current_x;
2181
2182 return SCIPnlpiOracleEvalConstraintValues(scip, nlpiproblem->oracle, x, g) == SCIP_OKAY;
2183}
2184
2185/** Method to return:
2186 * 1) The structure of the jacobian (if "values" is NULL)
2187 * 2) The values of the jacobian (if "values" is not NULL)
2188 */ /*lint -e{715}*/
2189bool ScipNLP::eval_jac_g(
2190 Index n, /**< number of variables */
2191 const Number* x, /**< point to evaluate */
2192 bool new_x, /**< whether some function evaluation method has been called for this point before */
2193 Index m, /**< number of constraints */
2194 Index nele_jac, /**< number of nonzero entries in jacobian */
2195 Index* iRow, /**< buffer to store row indices of nonzero jacobian entries, or NULL if values are requested */
2196 Index* jCol, /**< buffer to store column indices of nonzero jacobian entries, or NULL if values are requested */
2197 Number* values /**< buffer to store values of nonzero jacobian entries, or NULL if structure is requested */
2198 )
2199{ /*lint --e{715}*/
2200 assert(nlpiproblem != NULL);
2201 assert(nlpiproblem->oracle != NULL);
2202
2203 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2204 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2205
2206 if( values == NULL )
2207 { /* Ipopt wants to know sparsity structure */
2208 const int* jacoffset;
2209 const int* jaccol;
2210 int j;
2211 int i;
2212
2213 assert(iRow != NULL);
2214 assert(jCol != NULL);
2215
2216 if( SCIPnlpiOracleGetJacobianSparsity(scip, nlpiproblem->oracle, &jacoffset, &jaccol) != SCIP_OKAY )
2217 return false;
2218
2219 assert(jacoffset[0] == 0);
2220 assert(jacoffset[m] == nele_jac);
2221 j = jacoffset[0];
2222 for( i = 0; i < m; ++i )
2223 for( ; j < jacoffset[i+1]; ++j )
2224 iRow[j] = i;
2225
2227 }
2228 else
2229 {
2230 if( new_x )
2231 ++current_x;
2232 else
2233 {
2234 // pass new_x = TRUE to Jacobian eval iff we have not evaluated the constraint functions at this point yet
2235 new_x = last_g_eval_x < current_x;
2236 }
2237 // if we evaluate the Jacobian with new_x = true, then this will also evaluate the constraint functions
2238 // (and if we do with new_x = false, then we already have last_g_eval_x == current_x anyway)
2239 last_f_eval_x = current_x;
2240
2241 if( SCIPnlpiOracleEvalJacobian(scip, nlpiproblem->oracle, x, new_x, NULL, values) != SCIP_OKAY )
2242 return false;
2243 }
2244
2245 return true;
2246}
2247
2248/** Method to return:
2249 * 1) The structure of the hessian of the lagrangian (if "values" is NULL)
2250 * 2) The values of the hessian of the lagrangian (if "values" is not NULL)
2251 */ /*lint -e{715}*/
2252bool ScipNLP::eval_h(
2253 Index n, /**< number of variables */
2254 const Number* x, /**< point to evaluate */
2255 bool new_x, /**< whether some function evaluation method has been called for this point before */
2256 Number obj_factor, /**< weight for objective function */
2257 Index m, /**< number of constraints */
2258 const Number* lambda, /**< weights for constraint functions */
2259 bool new_lambda, /**< whether the hessian has been evaluated for these values of lambda before */
2260 Index nele_hess, /**< number of nonzero entries in hessian */
2261 Index* iRow, /**< buffer to store row indices of nonzero hessian entries, or NULL if values are requested */
2262 Index* jCol, /**< buffer to store column indices of nonzero hessian entries, or NULL if values are requested */
2263 Number* values /**< buffer to store values of nonzero hessian entries, or NULL if structure is requested */
2264 )
2265{ /*lint --e{715}*/
2266 assert(nlpiproblem != NULL);
2267 assert(nlpiproblem->oracle != NULL);
2268
2269 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2270 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2271
2272 if( values == NULL )
2273 { /* Ipopt wants to know sparsity structure */
2274 const int* heslagoffset;
2275 const int* heslagcol;
2276 int j;
2277 int i;
2278
2279 assert(iRow != NULL);
2280 assert(jCol != NULL);
2281
2283 return false;
2284
2285 assert(heslagoffset[0] == 0);
2287 j = heslagoffset[0];
2288 for( i = 0; i < n; ++i )
2289 for( ; j < heslagoffset[i+1]; ++j )
2290 iRow[j] = i;
2291
2293 }
2294 else
2295 {
2296 bool new_x_obj = new_x;
2297 bool new_x_cons = new_x;
2298 if( new_x )
2299 ++current_x;
2300 else
2301 {
2302 // pass new_x_obj = TRUE iff we have not evaluated the objective function at this point yet
2303 // pass new_x_cons = TRUE iff we have not evaluated the constraint functions at this point yet
2304 new_x_obj = last_f_eval_x < current_x;
2305 new_x_cons = last_g_eval_x < current_x;
2306 }
2307 // evaluating Hessians with new_x will also evaluate the functions itself
2308 last_f_eval_x = current_x;
2309 last_g_eval_x = current_x;
2310
2311 if( SCIPnlpiOracleEvalHessianLag(scip, nlpiproblem->oracle, x, new_x_obj, new_x_cons, obj_factor, lambda, values) != SCIP_OKAY )
2312 return false;
2313 }
2314
2315 return true;
2316}
2317
2318/** Method called by the solver at each iteration.
2319 *
2320 * Checks whether SCIP solve is interrupted, objlimit is reached, or fastfail is triggered.
2321 * Sets solution and termination status accordingly.
2322 */ /*lint -e{715}*/
2323bool ScipNLP::intermediate_callback(
2324 AlgorithmMode mode, /**< current mode of algorithm */
2325 Index iter, /**< current iteration number */
2326 Number obj_value, /**< current objective value */
2327 Number inf_pr, /**< current primal infeasibility */
2328 Number inf_du, /**< current dual infeasibility */
2329 Number mu, /**< current barrier parameter */
2330 Number d_norm, /**< current gradient norm */
2331 Number regularization_size,/**< current size of regularization */
2332 Number alpha_du, /**< current dual alpha */
2333 Number alpha_pr, /**< current primal alpha */
2334 Index ls_trials, /**< current number of linesearch trials */
2335 const IpoptData* ip_data, /**< pointer to Ipopt Data */
2336 IpoptCalculatedQuantities* ip_cq /**< pointer to current calculated quantities */
2337 )
2338{ /*lint --e{715}*/
2340 {
2341 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2342 nlpiproblem->termstat = SCIP_NLPTERMSTAT_INTERRUPT;
2343 return false;
2344 }
2345
2346 /* feasible point with objective value below lower objective limit -> stop */
2347 if( obj_value <= param.lobjlimit && inf_pr <= param.feastol )
2348 {
2349 nlpiproblem->solstat = SCIP_NLPSOLSTAT_FEASIBLE;
2350 nlpiproblem->termstat = SCIP_NLPTERMSTAT_LOBJLIMIT;
2351 return false;
2352 }
2353
2354 /* do convergence test if fastfail is enabled */
2355 if( param.fastfail >= SCIP_NLPPARAM_FASTFAIL_AGGRESSIVE )
2356 {
2357 int i;
2358
2359 if( iter == 0 )
2360 {
2361 conv_lastrestoiter = -1;
2362 }
2363 else if( mode == RestorationPhaseMode )
2364 {
2365 conv_lastrestoiter = iter;
2366 }
2367 else if( conv_lastrestoiter == iter-1 )
2368 {
2369 /* just switched back from restoration mode, reset dual reduction targets */
2370 for( i = 0; i < convcheck_nchecks; ++i )
2371 conv_dutarget[i] = convcheck_minred[i] * inf_du;
2372 }
2373
2374 if( iter == convcheck_startiter )
2375 {
2376 /* define initial targets and iteration limits */
2377 for( i = 0; i < convcheck_nchecks; ++i )
2378 {
2379 conv_prtarget[i] = convcheck_minred[i] * inf_pr;
2380 conv_dutarget[i] = convcheck_minred[i] * inf_du;
2381 conv_iterlim[i] = iter + convcheck_maxiter[i];
2382 }
2383 }
2384 else if( iter > convcheck_startiter )
2385 {
2386 /* check if we should stop */
2387 for( i = 0; i < convcheck_nchecks; ++i )
2388 {
2389 if( inf_pr <= conv_prtarget[i] )
2390 {
2391 /* sufficient reduction w.r.t. primal infeasibility target
2392 * reset target w.r.t. current infeasibilities
2393 */
2394 conv_prtarget[i] = convcheck_minred[i] * inf_pr;
2395 conv_dutarget[i] = convcheck_minred[i] * inf_du;
2396 conv_iterlim[i] = iter + convcheck_maxiter[i];
2397 }
2398 else if( iter >= conv_iterlim[i] )
2399 {
2400 /* we hit a limit, should we really stop? */
2401 SCIPdebugMsg(scip, "convcheck %d: inf_pr = %e > target %e; inf_du = %e target %e: ",
2402 i, inf_pr, conv_prtarget[i], inf_du, conv_dutarget[i]);
2403 if( mode == RegularMode && iter <= conv_lastrestoiter + convcheck_startiter )
2404 {
2405 /* if we returned from feasibility restoration recently, we allow some more iterations,
2406 * because Ipopt may go for optimality for some iterations, at the costs of infeasibility
2407 */
2408 SCIPdebugPrintf("continue, because restoration phase only %d iters ago\n", iter - conv_lastrestoiter);
2409 }
2410 else if( mode == RegularMode && inf_du <= conv_dutarget[i] && iter < conv_iterlim[i] + convcheck_maxiter[i] )
2411 {
2412 /* if dual reduction is sufficient, we allow for twice the number of iterations to reach primal infeas reduction */
2413 SCIPdebugPrintf("continue, because dual infeas. red. sufficient and only %d iters above limit\n", iter - conv_iterlim[i]);
2414 }
2415 else
2416 {
2417 SCIPdebugPrintf("abort solve\n");
2418 if( inf_pr <= param.feastol )
2419 nlpiproblem->solstat = SCIP_NLPSOLSTAT_FEASIBLE;
2420 else
2421 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2422 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OKAY;
2423 return false;
2424 }
2425 }
2426 }
2427 }
2428 }
2429
2430 return true;
2431}
2432
2433/** This method is called when the algorithm is complete so the TNLP can store/write the solution. */ /*lint -e{715}*/
2434void ScipNLP::finalize_solution(
2435 SolverReturn status, /**< solve and solution status */
2436 Index n, /**< number of variables */
2437 const Number* x, /**< primal solution values */
2438 const Number* z_L, /**< dual values of variable lower bounds */
2439 const Number* z_U, /**< dual values of variable upper bounds */
2440 Index m, /**< number of constraints */
2441 const Number* g, /**< values of constraints */
2442 const Number* lambda, /**< dual values of constraints */
2443 Number obj_value, /**< objective function value */
2444 const IpoptData* data, /**< pointer to Ipopt Data */
2445 IpoptCalculatedQuantities* cq /**< pointer to calculated quantities */
2446 )
2447{ /*lint --e{715}*/
2448 assert(nlpiproblem != NULL);
2449 assert(nlpiproblem->oracle != NULL);
2450
2451 assert(n == SCIPnlpiOracleGetNVars(nlpiproblem->oracle));
2452 assert(m == SCIPnlpiOracleGetNConstraints(nlpiproblem->oracle));
2453
2454 bool check_feasibility = false; // whether we should check x for feasibility, if not NULL
2455 switch( status )
2456 {
2457 case SUCCESS:
2458 nlpiproblem->solstat = SCIP_NLPSOLSTAT_LOCOPT;
2459 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OKAY;
2460 assert(x != NULL);
2461 break;
2462
2464 /* if stop at acceptable point, then dual infeasibility can be arbitrary large, so claim only feasibility */
2466 nlpiproblem->solstat = SCIP_NLPSOLSTAT_FEASIBLE;
2467 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OKAY;
2468 assert(x != NULL);
2469 break;
2470
2471 case MAXITER_EXCEEDED:
2472 check_feasibility = true;
2473 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2474 nlpiproblem->termstat = SCIP_NLPTERMSTAT_ITERLIMIT;
2475 break;
2476
2477 case CPUTIME_EXCEEDED:
2478#if IPOPT_VERSION_MAJOR > 3 || IPOPT_VERSION_MINOR >= 14
2479 case WALLTIME_EXCEEDED:
2480#endif
2481 check_feasibility = true;
2482 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2483 nlpiproblem->termstat = SCIP_NLPTERMSTAT_TIMELIMIT;
2484 break;
2485
2486 case STOP_AT_TINY_STEP:
2489 check_feasibility = true;
2490 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2491 nlpiproblem->termstat = SCIP_NLPTERMSTAT_NUMERICERROR;
2492 break;
2493
2495 nlpiproblem->solstat = SCIP_NLPSOLSTAT_LOCINFEASIBLE;
2496 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OKAY;
2497 break;
2498
2500 // status codes already set in intermediate_callback
2501 break;
2502
2503 case DIVERGING_ITERATES:
2504 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNBOUNDED;
2505 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OKAY;
2506 break;
2507
2508 // for the following status codes, if we get called here at all,
2509 // then Ipopt passes zeros for duals and activities!
2510 // (see https://github.com/coin-or/Ipopt/blob/stable/3.14/src/Interfaces/IpIpoptApplication.cpp#L885-L934)
2511
2513 // we can get this, if functions can still be evaluated, but are not differentiable
2514 // (so Ipopt couldn't check local optimality)
2515 // so we enable the check below for whether the point is feasible
2516 check_feasibility = true;
2517 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2518 nlpiproblem->termstat = SCIP_NLPTERMSTAT_EVALERROR;
2519 break;
2520
2522 case INTERNAL_ERROR:
2523 case INVALID_OPTION:
2524 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2525 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OTHER;
2526 break;
2527
2528 case OUT_OF_MEMORY:
2529 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2530 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OUTOFMEMORY;
2531 break;
2532
2533 default:
2534 SCIPerrorMessage("Ipopt returned with unknown solution status %d\n", status);
2535 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2536 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OTHER;
2537 break;
2538 }
2539
2540 assert(x != NULL);
2541 assert(lambda != NULL);
2542 assert(z_L != NULL);
2543 assert(z_U != NULL);
2544
2545 assert(nlpiproblem->solprimals != NULL);
2546
2547 if( nlpiproblem->soldualcons == NULL )
2548 {
2549 (void) SCIPallocBlockMemoryArray(scip, &nlpiproblem->soldualcons, m);
2550 }
2551 if( nlpiproblem->soldualvarlb == NULL )
2552 {
2553 (void) SCIPallocBlockMemoryArray(scip, &nlpiproblem->soldualvarlb, n);
2554 }
2555 if( nlpiproblem->soldualvarub == NULL )
2556 {
2557 (void) SCIPallocBlockMemoryArray(scip, &nlpiproblem->soldualvarub, n);
2558 }
2559 if( nlpiproblem->soldualcons == NULL || nlpiproblem->soldualvarlb == NULL || nlpiproblem->soldualvarub == NULL )
2560 {
2561 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2562 nlpiproblem->termstat = SCIP_NLPTERMSTAT_OUTOFMEMORY;
2563 return;
2564 }
2565
2566 BMScopyMemoryArray(nlpiproblem->solprimals, x, n);
2567 BMScopyMemoryArray(nlpiproblem->soldualcons, lambda, m);
2568 BMScopyMemoryArray(nlpiproblem->soldualvarlb, z_L, n);
2569 BMScopyMemoryArray(nlpiproblem->soldualvarub, z_U, n);
2570 nlpiproblem->solobjval = obj_value;
2571 nlpiproblem->solprimalvalid = true;
2572 nlpiproblem->solprimalgiven = false;
2573 nlpiproblem->soldualvalid = true;
2574 nlpiproblem->soldualgiven = false;
2575
2576 // get violations, there could be an evaluation error when doing so
2577#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2578 nlpiproblem->solboundviol = 0.0; // old Ipopt does not calculate bound violations, but for what it's worth, we have set bound_relax_factor=0 then
2579 if( cq == NULL )
2580 {
2581 // with old Ipopt, finalize_solution may be called with cq == NULL if all variables are fixed; we just skip the rest then
2582 nlpiproblem->solconsviol = 0.0;
2583 return;
2584 }
2585#else
2586 assert(cq != NULL);
2587 nlpiproblem->solboundviol = cq->unscaled_curr_orig_bounds_violation(Ipopt::NORM_MAX);
2588#endif
2589 try
2590 {
2591 nlpiproblem->solconsviol = cq->unscaled_curr_nlp_constraint_violation(Ipopt::NORM_MAX);
2592
2593 if( check_feasibility )
2594 {
2595 // we assume that check_feasibility has not been enabled if Ipopt claimed infeasibility, since we should not change solstatus to unknown then
2596 assert(nlpiproblem->solstat != SCIP_NLPSOLSTAT_LOCINFEASIBLE);
2597 if( MAX(nlpiproblem->solconsviol, nlpiproblem->solboundviol) <= param.feastol )
2598 nlpiproblem->solstat = SCIP_NLPSOLSTAT_FEASIBLE;
2599 else
2600 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2601 }
2602 }
2603 catch( const IpoptNLP::Eval_Error& exc )
2604 {
2605 SCIPdebugMsg(scip, "Eval error when checking constraint viol: %s\n", exc.Message().c_str());
2607 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2608 nlpiproblem->solconsviol = SCIP_INVALID;
2609 }
2610 catch(...)
2611 {
2612 /* with clang++, an IpoptNLP::Eval_Error wasn't catched by the catch-block above
2613 * I don't know why, but this should work around it
2614 */
2615 SCIPdebugMsg(scip, "Unknown exception when checking constraint viol\n");
2617 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2618 nlpiproblem->solconsviol = SCIP_INVALID;
2619 }
2620
2621 if( nlpiproblem->solstat == SCIP_NLPSOLSTAT_LOCINFEASIBLE )
2622 {
2623 assert(lambda != NULL);
2624 SCIP_Real tol;
2625 (void) nlpiproblem->ipopt->Options()->GetNumericValue("tol", tol, "");
2626
2627 // Jakobs paper ZR_20-20 says we should have lambda*g(x) + mu*h(x) > 0
2628 // if the NLP is min f(x) s.t. g(x) <= 0, h(x) = 0
2629 // we check this here and change solution status to unknown if the test fails
2630 bool infreasonable = true;
2631 SCIP_Real infproof = 0.0;
2632 for( int i = 0; i < m && infreasonable; ++i )
2633 {
2634 if( fabs(lambda[i]) < tol )
2635 continue;
2636 SCIP_Real side;
2637 if( lambda[i] < 0.0 )
2638 {
2639 // lhs <= g(x) should be active
2640 // in the NLP above, this should be lhs - g(x) <= 0 with negated dual
2641 // so this contributes -lambda*(lhs-g(x)) = lambda*(g(x)-side)
2642 side = SCIPnlpiOracleGetConstraintLhs(nlpiproblem->oracle, i);
2643 if( SCIPisInfinity(scip, -side) )
2644 {
2645 SCIPdebugMessage("inconsistent dual, lambda = %g, but lhs = %g\n", lambda[i], side);
2646 infreasonable = false;
2647 }
2648 }
2649 else
2650 {
2651 // g(x) <= rhs should be active
2652 // in the NLP above, this should be g(x) - rhs <= 0
2653 // so this contributes lambda*(g(x)-rhs)
2654 side = SCIPnlpiOracleGetConstraintRhs(nlpiproblem->oracle, i);
2655 if( SCIPisInfinity(scip, side) )
2656 {
2657 SCIPdebugMessage("inconsistent dual, lambda = %g, but rhs = %g\n", lambda[i], side);
2658 infreasonable = false;
2659 }
2660 }
2661
2662 // g(x) <= 0
2663 infproof += lambda[i] * (g[i] - side);
2664 // SCIPdebugMessage("cons %d lambda %g, slack %g\n", i, lambda[i], g[i] - side);
2665 }
2666 if( infreasonable )
2667 {
2668 SCIPdebugMessage("infproof = %g should be positive to be valid\n", infproof);
2669 if( infproof <= 0.0 )
2670 infreasonable = false;
2671 }
2672
2673 if( !infreasonable )
2674 {
2675 // change status to say we don't know
2676 nlpiproblem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
2677 }
2678 }
2679}
2680
2681/** Calls Lapacks Dsyev routine to compute eigenvalues and eigenvectors of a dense matrix.
2682 *
2683 * It's here, because we use Ipopt's C interface to Lapack.
2684 */
2686 SCIP_Bool computeeigenvectors,/**< should also eigenvectors should be computed ? */
2687 int N, /**< dimension */
2688 SCIP_Real* a, /**< matrix data on input (size N*N); eigenvectors on output if computeeigenvectors == TRUE */
2689 SCIP_Real* w /**< buffer to store eigenvalues (size N) */
2690 )
2691{
2692 int info;
2693
2694#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2695 IpLapackDsyev((bool)computeeigenvectors, N, a, N, w, info);
2696#else
2697 IpLapackSyev((bool)computeeigenvectors, N, a, N, w, info);
2698#endif
2699
2700 if( info != 0 )
2701 {
2702 SCIPerrorMessage("There was an error when calling DSYEV. INFO = %d\n", info);
2703 return SCIP_ERROR;
2704 }
2705
2706 return SCIP_OKAY;
2707}
2708
2709/** solves a linear problem of the form Ax = b for a regular matrix 3*3 A */
2710static
2712 SCIP_Real* A, /**< matrix data on input (size 3*3); filled column-wise */
2713 SCIP_Real* b, /**< right hand side vector (size 3) */
2714 SCIP_Real* x, /**< buffer to store solution (size 3) */
2715 SCIP_Bool* success /**< pointer to store if the solving routine was successful */
2716 )
2717{
2718 SCIP_Real Acopy[9];
2719 SCIP_Real bcopy[3];
2720 int pivotcopy[3];
2721 const int N = 3;
2722 int info;
2723
2724 assert(A != NULL);
2725 assert(b != NULL);
2726 assert(x != NULL);
2727 assert(success != NULL);
2728
2731
2732 /* compute the LU factorization */
2733#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2734 IpLapackDgetrf(N, Acopy, pivotcopy, N, info);
2735#else
2736 IpLapackGetrf(N, Acopy, pivotcopy, N, info);
2737#endif
2738
2739 if( info != 0 )
2740 {
2741 SCIPdebugMessage("There was an error when calling Dgetrf. INFO = %d\n", info);
2742 *success = FALSE;
2743 }
2744 else
2745 {
2746 *success = TRUE;
2747
2748 /* solve linear problem */
2749#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2751#else
2753#endif
2754
2755 /* copy the solution */
2757 }
2758
2759 return SCIP_OKAY;
2760}
2761
2762/** solves a linear problem of the form Ax = b for a regular matrix A
2763 *
2764 * Calls Lapacks DGETRF routine to calculate a LU factorization and uses this factorization to solve
2765 * the linear problem Ax = b.
2766 *
2767 * It's here, because we use Ipopt's C interface to Lapack.
2768 */
2770 int N, /**< dimension */
2771 SCIP_Real* A, /**< matrix data on input (size N*N); filled column-wise */
2772 SCIP_Real* b, /**< right hand side vector (size N) */
2773 SCIP_Real* x, /**< buffer to store solution (size N) */
2774 SCIP_Bool* success /**< pointer to store if the solving routine was successful */
2775 )
2776{
2777 SCIP_Real* Acopy;
2778 SCIP_Real* bcopy;
2779 int* pivotcopy;
2780 int info;
2781
2782 assert(N > 0);
2783 assert(A != NULL);
2784 assert(b != NULL);
2785 assert(x != NULL);
2786 assert(success != NULL);
2787
2788 /* call solveLinearProb3() for performance reasons */
2789 if( N == 3 )
2790 {
2792 return SCIP_OKAY;
2793 }
2794
2795 Acopy = NULL;
2796 bcopy = NULL;
2797 pivotcopy = NULL;
2798
2802
2803 /* compute the LU factorization */
2804#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2805 IpLapackDgetrf(N, Acopy, pivotcopy, N, info);
2806#else
2807 IpLapackGetrf(N, Acopy, pivotcopy, N, info);
2808#endif
2809
2810 if( info != 0 )
2811 {
2812 SCIPdebugMessage("There was an error when calling Dgetrf. INFO = %d\n", info);
2813 *success = FALSE;
2814 }
2815 else
2816 {
2817 *success = TRUE;
2818
2819 /* solve linear problem */
2820#if IPOPT_VERSION_MAJOR == 3 && IPOPT_VERSION_MINOR < 14
2822#else
2824#endif
2825
2826 /* copy the solution */
2828 }
2829
2833
2834 return SCIP_OKAY;
2835}
SCIP_VAR * w
SCIP_VAR * a
SCIP_VAR ** b
SCIP_VAR ** x
#define SCIP_MAXSTRLEN
Definition def.h:302
#define SCIP_INVALID
Definition def.h:206
#define SCIP_ALLOC(x)
Definition def.h:399
#define SCIP_Real
Definition def.h:186
#define TRUE
Definition def.h:95
#define FALSE
Definition def.h:96
#define SCIP_CALL_ABORT(x)
Definition def.h:367
#define SCIP_CALL(x)
Definition def.h:388
methods to interpret (evaluate) an expression "fast"
int SCIPgetSubscipDepth(SCIP *scip)
Definition scip_copy.c:2605
SCIP_EXPRINTCAPABILITY SCIPexprintGetCapability(void)
void SCIPinfoMessage(SCIP *scip, FILE *file, const char *formatstr,...)
SCIP_MESSAGEHDLR * SCIPgetMessagehdlr(SCIP *scip)
#define SCIPdebugMsg
SCIP_RETCODE SCIPincludeNlpSolverIpopt(SCIP *scip)
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveValue(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Real *objval)
SCIP_RETCODE SCIPnlpiOracleChgLinearCoefs(SCIP *scip, SCIP_NLPIORACLE *oracle, int considx, int nentries, const int *varidxs, const SCIP_Real *newcoefs)
SCIP_RETCODE SCIPnlpiOracleChgVarBounds(SCIP *scip, SCIP_NLPIORACLE *oracle, int nvars, const int *indices, const SCIP_Real *lbs, const SCIP_Real *ubs)
SCIP_RETCODE SCIPnlpiOracleAddConstraints(SCIP *scip, SCIP_NLPIORACLE *oracle, int nconss, const SCIP_Real *lhss, const SCIP_Real *rhss, const int *nlininds, int *const *lininds, SCIP_Real *const *linvals, SCIP_EXPR **exprs, const char **consnames)
SCIP_Bool SCIPnlpiOracleIsConstraintNonlinear(SCIP_NLPIORACLE *oracle, int considx)
SCIP_RETCODE SCIPnlpiOracleDelVarSet(SCIP *scip, SCIP_NLPIORACLE *oracle, int *delstats)
SCIP_RETCODE SCIPnlpiOracleEvalConstraintValues(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Real *convals)
SCIP_RETCODE SCIPnlpiOracleCreate(SCIP *scip, SCIP_NLPIORACLE **oracle)
Definition nlpioracle.c:983
SCIP_RETCODE SCIPnlpiOracleGetJacobianSparsity(SCIP *scip, SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
void SCIPnlpiOracleGetVarCounts(SCIP *scip, SCIP_NLPIORACLE *oracle, const int **lincounts, const int **nlcounts)
SCIP_RETCODE SCIPnlpiOracleGetHessianLagSparsity(SCIP *scip, SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
char * SCIPnlpiOracleGetConstraintName(SCIP_NLPIORACLE *oracle, int considx)
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveGradient(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *objval, SCIP_Real *objgrad)
SCIP_RETCODE SCIPnlpiOracleResetEvalTime(SCIP *scip, SCIP_NLPIORACLE *oracle)
SCIP_RETCODE SCIPnlpiOracleSetObjective(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real constant, int nlin, const int *lininds, const SCIP_Real *linvals, SCIP_EXPR *expr)
SCIP_Real SCIPnlpiOracleGetConstraintRhs(SCIP_NLPIORACLE *oracle, int considx)
SCIP_Real SCIPnlpiOracleGetEvalTime(SCIP *scip, SCIP_NLPIORACLE *oracle)
SCIP_RETCODE SCIPnlpiOracleChgConsSides(SCIP *scip, SCIP_NLPIORACLE *oracle, int nconss, const int *indices, const SCIP_Real *lhss, const SCIP_Real *rhss)
SCIP_Real SCIPnlpiOracleGetConstraintLhs(SCIP_NLPIORACLE *oracle, int considx)
SCIP_RETCODE SCIPnlpiOracleAddVars(SCIP *scip, SCIP_NLPIORACLE *oracle, int nvars, const SCIP_Real *lbs, const SCIP_Real *ubs, const char **varnames)
int SCIPnlpiOracleGetNVars(SCIP_NLPIORACLE *oracle)
int SCIPnlpiOracleGetNConstraints(SCIP_NLPIORACLE *oracle)
SCIP_EXPRINTCAPABILITY SCIPnlpiOracleGetEvalCapability(SCIP *scip, SCIP_NLPIORACLE *oracle)
SCIP_Real SCIPnlpiOracleGetObjectiveConstant(SCIP_NLPIORACLE *oracle)
SCIP_RETCODE SCIPnlpiOracleEvalHessianLag(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx_obj, SCIP_Bool isnewx_cons, SCIP_Real objfactor, const SCIP_Real *lambda, SCIP_Real *hessian)
SCIP_Bool SCIPnlpiOracleIsVarNonlinear(SCIP *scip, SCIP_NLPIORACLE *oracle, int varidx)
SCIP_RETCODE SCIPnlpiOracleEvalJacobian(SCIP *scip, SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *convals, SCIP_Real *jacobi)
SCIP_RETCODE SCIPnlpiOracleDelConsSet(SCIP *scip, SCIP_NLPIORACLE *oracle, int *delstats)
SCIP_RETCODE SCIPnlpiOracleSetProblemName(SCIP *scip, SCIP_NLPIORACLE *oracle, const char *name)
SCIP_RETCODE SCIPnlpiOracleChgObjConstant(SCIP *scip, SCIP_NLPIORACLE *oracle, SCIP_Real objconstant)
char ** SCIPnlpiOracleGetVarNames(SCIP_NLPIORACLE *oracle)
const SCIP_Real * SCIPnlpiOracleGetVarLbs(SCIP_NLPIORACLE *oracle)
const SCIP_Real * SCIPnlpiOracleGetVarUbs(SCIP_NLPIORACLE *oracle)
SCIP_RETCODE SCIPnlpiOracleFree(SCIP *scip, SCIP_NLPIORACLE **oracle)
const char * SCIPnlpiOracleGetProblemName(SCIP_NLPIORACLE *oracle)
SCIP_RETCODE SCIPnlpiOracleChgExpr(SCIP *scip, SCIP_NLPIORACLE *oracle, int considx, SCIP_EXPR *expr)
const char * SCIPgetSolverNameIpopt(void)
SCIP_RETCODE SCIPcallLapackDsyevIpopt(SCIP_Bool computeeigenvectors, int N, SCIP_Real *a, SCIP_Real *w)
SCIP_RETCODE SCIPsolveLinearEquationsIpopt(int N, SCIP_Real *A, SCIP_Real *b, SCIP_Real *x, SCIP_Bool *success)
void * SCIPgetNlpiOracleIpopt(SCIP_NLPIPROBLEM *nlpiproblem)
SCIP_Bool SCIPisIpoptAvailableIpopt(void)
const char * SCIPgetSolverDescIpopt(void)
SCIP_RETCODE SCIPaddIntParam(SCIP *scip, const char *name, const char *desc, int *valueptr, SCIP_Bool isadvanced, int defaultvalue, int minvalue, int maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition scip_param.c:83
SCIP_PARAM * SCIPgetParam(SCIP *scip, const char *name)
Definition scip_param.c:234
SCIP_RETCODE SCIPaddStringParam(SCIP *scip, const char *name, const char *desc, char **valueptr, SCIP_Bool isadvanced, const char *defaultvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition scip_param.c:194
SCIP_RETCODE SCIPaddRealParam(SCIP *scip, const char *name, const char *desc, SCIP_Real *valueptr, SCIP_Bool isadvanced, SCIP_Real defaultvalue, SCIP_Real minvalue, SCIP_Real maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition scip_param.c:139
SCIP_RETCODE SCIPincludeExternalCodeInformation(SCIP *scip, const char *name, const char *description)
#define SCIPallocBlockMemoryArray(scip, ptr, num)
Definition scip_mem.h:93
#define SCIPreallocBlockMemoryArray(scip, ptr, oldnum, newnum)
Definition scip_mem.h:99
#define SCIPfreeBlockMemoryArrayNull(scip, ptr, num)
Definition scip_mem.h:111
SCIP_RETCODE SCIPincludeNlpi(SCIP *scip, 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_NLPICHGEXPR((*nlpichgexpr)), 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_NLPIDATA *nlpidata)
Definition scip_nlpi.c:107
SCIP_NLPIDATA * SCIPnlpiGetData(SCIP_NLPI *nlpi)
Definition nlpi.c:712
SCIP_Bool SCIPisSolveInterrupted(SCIP *scip)
SCIP_Real SCIPinfinity(SCIP *scip)
SCIP_Real SCIPfloor(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisInfinity(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Real SCIPrandomGetReal(SCIP_RANDNUMGEN *randnumgen, SCIP_Real minrandval, SCIP_Real maxrandval)
Definition misc.c:10041
int SCIPsnprintf(char *t, int len, const char *s,...)
Definition misc.c:10788
return SCIP_OKAY
SCIPfreeRandom(scip, &heurdata->randnumgen)
SCIP_Real objval
SCIPcreateRandom(scip, &heurdata->randnumgen, DEFAULT_RANDSEED, TRUE))
assert(minobj< SCIPgetCutoffbound(scip))
int nvars
static const char * paramname[]
Definition lpi_msk.c:5096
#define NULL
Definition lpi_spx1.cpp:161
#define BMSduplicateMemoryArray(ptr, source, num)
Definition memory.h:145
#define BMSallocMemoryArray(ptr, num)
Definition memory.h:125
#define BMSfreeMemoryArray(ptr)
Definition memory.h:149
#define BMScopyMemoryArray(ptr, source, num)
Definition memory.h:136
void SCIPmessageVPrintError(const char *formatstr, va_list ap)
Definition message.c:804
void SCIPmessageVPrintInfo(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr, va_list ap)
Definition message.c:608
void SCIPmessagePrintErrorHeader(const char *sourcefile, int sourceline)
Definition message.c:777
#define NLPI_PRIORITY
#define MAXPERTURB
static const SCIP_Real convcheck_minred[convcheck_nchecks]
static const char * ipopt_int_params[]
integer parameters of Ipopt to make available via SCIP parameters
#define NLPI_NAME
static const int convcheck_startiter
static const int convcheck_nchecks
#define DEFAULT_RANDSEED
static SCIP_RETCODE solveLinearProb3(SCIP_Real *A, SCIP_Real *b, SCIP_Real *x, SCIP_Bool *success)
#define FEASTOLFACTOR
static SCIP_RETCODE handleNlpParam(SCIP *scip, SCIP_NLPIDATA *nlpidata, SCIP_NLPIPROBLEM *nlpiproblem, const SCIP_NLPPARAM param)
pass NLP solve parameters to Ipopt
static const char * ipopt_string_params[]
string parameters of Ipopt to make available via SCIP parameters
static SCIP_RETCODE ensureStartingPoint(SCIP *scip, SCIP_NLPIPROBLEM *problem, SCIP_Bool &warmstart)
static const int convcheck_maxiter[convcheck_nchecks]
static void invalidateSolved(SCIP_NLPIPROBLEM *problem)
#define NLPI_DESC
static void invalidateSolution(SCIP_NLPIPROBLEM *problem)
Ipopt NLP interface.
methods to store an NLP and request function, gradient, and Hessian values
char * SCIPparamGetString(SCIP_PARAM *param)
Definition paramset.c:911
int SCIPparamGetInt(SCIP_PARAM *param)
Definition paramset.c:734
int SCIPparamGetIntDefault(SCIP_PARAM *param)
Definition paramset.c:770
#define SCIPerrorMessage
Definition pub_message.h:64
#define SCIPdebugMessage
Definition pub_message.h:96
#define SCIPdebugPrintf
Definition pub_message.h:99
public data structures and miscellaneous methods
public methods for handling parameter settings
public methods for problem copies
general public methods
public methods for memory management
public methods for message handling
public methods for NLPI solver interfaces
public methods for numerical tolerances
public methods for SCIP parameter handling
public methods for random numbers
public solving methods
SCIP_Real timelimit
Definition type_nlpi.h:72
SCIP_Real feastol
Definition type_nlpi.h:69
SCIP_Bool expectinfeas
Definition type_nlpi.h:76
SCIP_Bool warmstart
Definition type_nlpi.h:77
SCIP_Real opttol
Definition type_nlpi.h:70
SCIP_Real solvertol
Definition type_nlpi.h:71
SCIP_NLPPARAM_FASTFAIL fastfail
Definition type_nlpi.h:75
unsigned short verblevel
Definition type_nlpi.h:74
SCIP_Real solconsviol
SCIP_Real * soldualcons
SCIP_Real solboundviol
SmartPtr< IpoptApplication > ipopt
SCIP_NLPIORACLE * oracle
SCIP_RANDNUMGEN * randnumgen
SCIP_Real * soldualvarlb
SCIP_NLPTERMSTAT termstat
SmartPtr< ScipNLP > nlp
SCIP_Real solobjval
SCIP_NLPSOLSTAT solstat
SCIP_Real * soldualvarub
SCIP_Real * solprimals
#define MAX(x, y)
Definition tclique_def.h:92
#define SCIP_EXPRINTCAPABILITY_GRADIENT
#define SCIP_EXPRINTCAPABILITY_HESSIAN
#define SCIP_EXPRINTCAPABILITY_FUNCVALUE
unsigned int SCIP_EXPRINTCAPABILITY
#define SCIP_DECL_NLPISOLVE(x)
Definition type_nlpi.h:497
#define SCIP_DECL_NLPICHGLINEARCOEFS(x)
Definition type_nlpi.h:432
#define SCIP_DECL_NLPICHGOBJCONSTANT(x)
Definition type_nlpi.h:463
#define SCIP_NLPPARAM_PRINT(param)
Definition type_nlpi.h:142
#define SCIP_DECL_NLPIGETSOLUTION(x)
Definition type_nlpi.h:545
#define SCIP_DECL_NLPISETOBJECTIVE(x)
Definition type_nlpi.h:344
#define SCIP_DECL_NLPICREATEPROBLEM(x)
Definition type_nlpi.h:255
#define SCIP_DECL_NLPIGETSTATISTICS(x)
Definition type_nlpi.h:562
#define SCIP_DECL_NLPIDELCONSSET(x)
Definition type_nlpi.h:415
#define SCIP_DECL_NLPICHGCONSSIDES(x)
Definition type_nlpi.h:383
#define SCIP_DECL_NLPIDELVARSET(x)
Definition type_nlpi.h:400
#define SCIP_DECL_NLPICHGEXPR(x)
Definition type_nlpi.h:449
#define SCIP_DECL_NLPIADDVARS(x)
Definition type_nlpi.h:297
@ SCIP_NLPPARAM_FASTFAIL_OFF
Definition type_nlpi.h:58
@ SCIP_NLPPARAM_FASTFAIL_AGGRESSIVE
Definition type_nlpi.h:60
enum SCIP_NlpSolStat SCIP_NLPSOLSTAT
Definition type_nlpi.h:168
#define SCIP_DECL_NLPISETINITIALGUESS(x)
Definition type_nlpi.h:481
#define SCIP_DECL_NLPIFREEPROBLEM(x)
Definition type_nlpi.h:267
@ SCIP_NLPTERMSTAT_OKAY
Definition type_nlpi.h:173
@ SCIP_NLPTERMSTAT_TIMELIMIT
Definition type_nlpi.h:174
@ SCIP_NLPTERMSTAT_NUMERICERROR
Definition type_nlpi.h:178
@ SCIP_NLPTERMSTAT_OTHER
Definition type_nlpi.h:182
@ SCIP_NLPTERMSTAT_EVALERROR
Definition type_nlpi.h:179
@ SCIP_NLPTERMSTAT_LOBJLIMIT
Definition type_nlpi.h:176
@ SCIP_NLPTERMSTAT_ITERLIMIT
Definition type_nlpi.h:175
@ SCIP_NLPTERMSTAT_OUTOFMEMORY
Definition type_nlpi.h:180
@ SCIP_NLPTERMSTAT_INTERRUPT
Definition type_nlpi.h:177
#define SCIP_DECL_NLPICOPY(x)
Definition type_nlpi.h:215
#define SCIP_DECL_NLPIGETSOLVERPOINTER(x)
Definition type_nlpi.h:243
#define SCIP_DECL_NLPIGETSOLSTAT(x)
Definition type_nlpi.h:511
#define SCIP_DECL_NLPICHGVARBOUNDS(x)
Definition type_nlpi.h:364
#define SCIP_DECL_NLPIGETPROBLEMPOINTER(x)
Definition type_nlpi.h:282
#define SCIP_DECL_NLPIFREE(x)
Definition type_nlpi.h:225
#define SCIP_DECL_NLPIADDCONSTRAINTS(x)
Definition type_nlpi.h:320
@ SCIP_NLPSOLSTAT_UNBOUNDED
Definition type_nlpi.h:165
@ SCIP_NLPSOLSTAT_LOCINFEASIBLE
Definition type_nlpi.h:163
@ SCIP_NLPSOLSTAT_FEASIBLE
Definition type_nlpi.h:162
@ SCIP_NLPSOLSTAT_LOCOPT
Definition type_nlpi.h:161
@ SCIP_NLPSOLSTAT_UNKNOWN
Definition type_nlpi.h:166
#define SCIP_DECL_NLPIGETTERMSTAT(x)
Definition type_nlpi.h:524
enum SCIP_NlpTermStat SCIP_NLPTERMSTAT
Definition type_nlpi.h:194
struct SCIP_NlpiData SCIP_NLPIDATA
Definition type_nlpi.h:52
@ SCIP_NOMEMORY
@ SCIP_ERROR
enum SCIP_Retcode SCIP_RETCODE