Scippy

SCIP

Solving Constraint Integer Programs

nlpi_filtersqp.c
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2017 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not email to scip@zib.de. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file nlpi_filtersqp.c
17  * @ingroup NLPIS
18  * @brief filterSQP NLP interface
19  * @author Stefan Vigerske
20  *
21  * @todo scaling
22  */
23 
24 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
25 
26 #include <string.h>
27 #if defined(_WIN32)
28 #include <windows.h>
29 #else
30 #include <sys/time.h>
31 #endif
32 
33 #ifndef NPARASCIP
34 #include <pthread.h>
35 #endif
36 
37 #include "scip/misc.h"
38 #include "nlpi/nlpi_filtersqp.h"
39 #include "nlpi/nlpi.h"
40 #include "nlpi/nlpioracle.h"
41 
42 #define NLPI_NAME "filtersqp" /* short concise name of solver */
43 #define NLPI_DESC "Sequential Quadratic Programming trust region solver by R. Fletcher and S. Leyffer" /* description of solver */
44 #define NLPI_PRIORITY -1000 /* priority of NLP solver */
45 
46 #define RANDSEED 26051979 /**< initial random seed */
47 #define MAXPERTURB 0.01 /**< maximal perturbation of bounds in starting point heuristic */
48 #define MAXNRUNS 3 /**< maximal number of FilterSQP calls per NLP solve (several calls if increasing workspace or decreasing eps) */
49 #define WORKSPACEGROWTHFACTOR 2 /**< factor by which to increase worksapce */
50 #define MINEPS 1e-14 /**< minimal FilterSQP epsilon */
51 #define OPTTOLFACTOR 0.5 /**< factor to apply to optimality tolerance, because FilterSQP do scaling */
52 #define DEFAULT_LOBJLIM (real)(-1e100) /**< default lower objective limit (should mean "unlimited") */
53 #define DEFAULT_FEASOPTTOL 1e-6 /**< default feasibility and optimality tolerance */
54 #define DEFAULT_MAXITER 3000 /**< default iteration limit */
55 
56 /*
57  * Data structures
58  */
59 
60 typedef int fint;
61 typedef double real;
62 typedef long ftnlen;
63 
64 struct SCIP_Time
65 {
66  time_t sec; /**< seconds part of time since epoch */
67  long usec; /**< micro-seconds part of time */
68 };
69 typedef struct SCIP_Time SCIP_TIME;
70 
71 struct SCIP_NlpiData
72 {
73  BMS_BLKMEM* blkmem; /**< block memory */
74  SCIP_MESSAGEHDLR* messagehdlr; /**< message handler */
75  SCIP_Real infinity; /**< initial value for infinity */
76  SCIP_RANDNUMGEN* randnumgen; /**< random number generator, if we have to make up a starting point */
77  SCIP_TIME starttime; /**< time at start of solve */
78 };
79 
80 struct SCIP_NlpiProblem
81 {
82  SCIP_NLPIORACLE* oracle; /**< Oracle-helper to store and evaluate NLP */
83 
84  int varssize; /**< length of variables-related arrays, if allocated */
85  int conssize; /**< length of constraints-related arrays, if allocated */
86 
87  SCIP_Real* initguess; /**< initial values for primal variables, or NULL if not known, size varssize */
88  SCIP_Bool warmstart; /**< whether we could warmstart the next solve */
89 
90  SCIP_Real* primalvalues; /**< primal values of variables in solution, size varssize */
91  SCIP_Real* consdualvalues; /**< dual values of constraints in solution, size conssize */
92  SCIP_Real* varlbdualvalues; /**< dual values of variable lower bounds in solution, size varssize */
93  SCIP_Real* varubdualvalues; /**< dual values of variable upper bounds in solution, size varssize */
94 
95  SCIP_NLPSOLSTAT solstat; /**< solution status from last NLP solve */
96  SCIP_NLPTERMSTAT termstat; /**< termination status from last NLP solve */
97  SCIP_Real solvetime; /**< time spend for last NLP solve */
98  int niterations; /**< number of iterations for last NLP solve */
99 
100  SCIP_Bool fromscratch; /**< value of fromscratch parameter */
101  fint istat[14]; /**< integer solution statistics from last FilterSQP call */
102  real rstat[7]; /**< real solution statistics from last FilterSQP call */
103  real feastol; /**< user-given feasibility tolerance */
104  real opttol; /**< user-given optimality tolerance */
105  real fmin; /**< lower bound on objective value */
106  SCIP_Real maxtime; /**< time limit */
107  fint maxiter; /**< iteration limit */
108  fint iprint; /**< print verbosity level */
109 
110  /* cached FilterSQP data */
111  real* x; /**< variable values, size varssize */
112  real* c; /**< constraint value, size conssize */
113  real* lam; /**< duals, size varssize + conssize */
114  real* bl; /**< variable lower bounds and constraint lhs, size varssize + conssize */
115  real* bu; /**< variable upper bounds and constraint rhs, size varssize + conssize */
116  real* s; /**< scaling factors, size varssize + conssize */
117  char* cstype; /**< constraint linearity, size conssize */
118  real* a; /**< gradients values, size la[0]-1 */
119  fint* la; /**< gradients indices, size lasize */
120  int lasize; /**< length of la array */
121  fint* hessiannz; /**< nonzero information about Hessian, size hessiannzsize */
122  int hessiannzsize;/**< length of hessiannz array */
123  real* ws; /**< real workspace, size mxwk */
124  fint* lws; /**< integer workspace, size mxiwk */
125  fint mxwk; /**< size of real workspace */
126  fint mxiwk; /**< size of integer workspace */
127  real* evalbuffer; /**< buffer to cache evaluation results before passing it to FilterSQP, size evalbufsize */
128  int evalbufsize; /**< size of evaluation buffer */
129 };
130 
131 /*
132  * Local methods
133  */
134 
135 #ifdef FNAME_LCASE_DECOR
136 # define F77_FUNC(name,NAME) name ## _
137 #endif
138 #ifdef FNAME_UCASE_DECOR
139 # define F77_FUNC(name,NAME) NAME ## _
140 #endif
141 #ifdef FNAME_LCASE_NODECOR
142 # define F77_FUNC(name,NAME) name
143 #endif
144 #ifdef FNAME_UCASE_NODECOR
145 # define F77_FUNC(name,NAME) NAME
146 #endif
147 
148 /** FilterSQP main routine.
149  *
150  * Array a has length nnza, which is the number of nonzeros in the gradient of the objective and the Jacobian.
151  * The first entries of a is the objective gradient, next are the gradients of the constraints.
152  *
153  * Array la has length lamax, which is at least nnza+m+2.
154  * la contains the index information of a row-oriented sparse matrix storage. It stores the number of nonzeros, the column indices, and the row starts:
155  * la[0] must be set to nnza+1.
156  * la[1]..la[nnza] are indices of the variables corresponding to the entries in a (colidx).
157  * la[nnza+1]..la[nnza+1+m] contain the index where each row starts in a and la (rowstart).
158  */
159 void F77_FUNC(filtersqp,FILTERSQP)(
160  fint* n, /**< number of variables */
161  fint* m, /**< number of constraints (excluding simple bounds) */
162  fint* kmax, /**< maximum size of null-space (at most n) */
163  fint* maxa, /**< maximum number of nonzero entries allowed in Jacobian */
164  fint* maxf, /**< maximum size of the filter - typically 100 */
165  fint* mlp, /**< maximum level parameter for resolving degeneracy in bqpd - typically 100 */
166  fint* mxwk, /**< maximum size of real workspace ws */
167  fint* mxiwk, /**< maximum size of integer workspace lws */
168  fint* iprint, /**< print flag: 0 is quiet, higher is more */
169  fint* nout, /**< output channel - 6 for screen */
170  fint* ifail, /**< fail flag and warmstart indicator */
171  real* rho, /**< initial trust-region radius - default 10 */
172  real* x, /**< starting point and final solution (array of length n) */
173  real* c, /**< final values of general constraints (array of length m) */
174  real* f, /**< final objective value */
175  real* fmin, /**< lower bound on objective value (as termination criteria) */
176  real* bl, /**< lower bounds of variables and constraints (array of length n+m) */
177  real* bu, /**< upper bounds of variables and constraints (array of length n+m) */
178  real* s, /**< scaling factors (array of length n+m) */
179  real* a, /**< objective gradient (always dense) and Jacobian (sparse or dense) entries */
180  fint* la, /**< column indices and length of rows of entries in a (if sparse) or leading dimension of a (if dense) */
181  real* ws, /**< real workspace (array of length mxwk) */
182  fint* lws, /**< integer workspace (array of length mxiwk) */
183  real* lam, /**< Lagrangian multipliers of simple bounds and general constraints (array of length n+m) */
184  char* cstype, /**< indicator whether constraint is linear ('L') or nonlinear ('N') (array of length m) */
185  real* user, /**< real workspace passed through to user routines */
186  fint* iuser, /**< integer workspace passed through to user routines */
187  fint* maxiter, /**< iteration limit for SQP solver */
188  fint* istat, /**< integer space for solution statistics (array of length 14) */
189  real* rstat, /**< real space for solution statitics (array of length 7) */
190  ftnlen cstype_len /**< 1 ??? */
191  );
192 
193 void F77_FUNC(objfun,OBJFUN)(real *x, fint *n, real *f, real *user, fint *iuser,
194  fint *errflag);
195 
196 void F77_FUNC(confun,CONFUN)(real *x, fint *n, fint *m, real *c, real *a, fint *la,
197  real *user, fint *iuser, fint *errflag);
198 
199 /* TODO what are the arguments of this function and does it need to be implemented?
200  * it's not in the filterSQP manual, but its an undefined symbol in the filterSQP library
201 void F77_FUNC(objgrad,OBJGRAD)(fint *, fint *, fint *, real *, real *, fint *, fint
202  *, real *, fint *, fint *);
203 */
204 void F77_FUNC(objgrad,OBJGRAD)(void);
205 
206 void F77_FUNC(gradient,GRADIENT)(fint *n, fint *m, fint *mxa, real *x, real *a, fint *la,
207  fint *maxa, real *user, fint *iuser, fint *errflag);
208 
209 void F77_FUNC(hessian,HESSIAN)(real *x, fint *n, fint *m, fint *phase, real *lam,
210  real *ws, fint *lws, real *user, fint *iuser,
211  fint *l_hess, fint *li_hess, fint *errflag);
212 
213 /** common block for problemname */
214 extern struct
215 {
216  fint char_l;
217  char pname[10];
218 } F77_FUNC(cpname,CPNAME);
219 /*lint -esym(752,cpname_) -esym(754,char_l) -esym(754,pname) */
220 
221 /** common block for Hessian storage set to 0, i.e. NO Hessian */
222 extern struct
223 {
225 } F77_FUNC(hessc,HESSC);
226 /*lint -esym(754,phr) -esym(754,phc) */
227 
228 /** common block for upper bound on filter */
229 extern struct
230 {
232 } F77_FUNC(ubdc,UBDC);
233 
234 /** common block for infinity & epsilon */
235 extern struct
236 {
238 } F77_FUNC(nlp_eps_inf,NLP_EPS_INF);
239 
240 /** common block for printing from QP solver */
241 extern struct
242 {
244 } F77_FUNC(bqpd_count,BQPD_COUNT);
245 /*lint -esym(752,bqpd_count_) -esym(754,n_bqpd_calls) -esym(754,n_bqpd_prfint) */
246 
247 /** common for scaling: scale_mode = 0 (none), 1 (variables), 2 (vars+cons) */
248 extern struct
249 {
251 } F77_FUNC(scalec,SCALEC);
252 /*lint -esym(754,phe) */
253 
254 #ifndef NPARASCIP
255 static pthread_mutex_t filtersqpmutex = PTHREAD_MUTEX_INITIALIZER;
256 #endif
257 
258 static
260 {
261  SCIP_TIME t;
262 #ifndef _WIN32
263  struct timeval tp; /*lint !e86*/
264 #endif
265 
266 #ifdef _WIN32
267  t.sec = time(NULL);
268  t.usec = 0;
269 
270 #else
271  (void) gettimeofday(&tp, NULL);
272  t.sec = tp.tv_sec;
273  t.usec = tp.tv_usec;
274 #endif
275 
276  return t;
277 }
278 
279 /* gives time since starttime (in problem) */
280 static
282  SCIP_NLPIDATA* nlpidata /**< NLPI data */
283  )
284 {
285  SCIP_TIME now;
286 
287  assert(nlpidata != NULL);
288 
289  now = gettime();
290 
291  /* now should be after startime */
292  assert(now.sec >= nlpidata->starttime.sec);
293  assert(now.sec > nlpidata->starttime.sec || now.usec >= nlpidata->starttime.usec);
294 
295  return (SCIP_Real)(now.sec - nlpidata->starttime.sec) + 1e-6 * (SCIP_Real)(now.usec - nlpidata->starttime.usec);
296 }
297 
298 static
300  SCIP_NLPIDATA* nlpidata, /**< NLPI data */
301  SCIP_NLPIPROBLEM* nlpiproblem /**< NLPI problem */
302  )
303 {
304  if( nlpiproblem->maxtime == DBL_MAX ) /*lint !e777 */
305  return FALSE;
306 
307  return timeelapsed(nlpidata) >= nlpiproblem->maxtime;
308 }
309 
310 /** Objective function evaluation */
311 void F77_FUNC(objfun,OBJFUN)(
312  real* x, /**< value of current variables (array of length n) */
313  fint* n, /**< number of variables */
314  real* f, /**< buffer to store value of objective function */
315  real* user, /**< user real workspace */
316  fint* iuser, /**< user integer workspace */
317  fint* errflag /**< set to 1 if arithmetic exception occurs, otherwise 0 */
318  )
319 { /*lint --e{715} */
320  SCIP_NLPIPROBLEM* problem;
321  real val;
322 
323  problem = (SCIP_NLPIPROBLEM*)(void*)iuser;
324  assert(problem != NULL);
325 
326  if( timelimitreached((SCIP_NLPIDATA*)(void*)user, problem) )
327  {
328  SCIPdebugMessage("timelimit reached, issuing arithmetic exception in objfun\n");
329  *errflag = 1;
330  return;
331  }
332 
333  *errflag = 1;
334  if( SCIPnlpiOracleEvalObjectiveValue(problem->oracle, x, &val) == SCIP_OKAY && SCIPisFinite(val) )
335  {
336  *errflag = 0;
337  *f = val;
338  }
339  else
340  {
341  SCIPdebugMessage("arithmetic exception in objfun\n");
342  }
343 }
344 
345 /** Constraint functions evaluation */
346 void F77_FUNC(confun,CONFUN)(
347  real* x, /**< value of current variables (array of length n) */
348  fint* n, /**< number of variables */
349  fint* m, /**< number of constraints */
350  real* c, /**< buffer to store values of constraints (array of length m) */
351  real* a, /**< Jacobian matrix entries */
352  fint* la, /**< Jacobian index information */
353  real* user, /**< user real workspace */
354  fint* iuser, /**< user integer workspace */
355  fint* errflag /**< set to 1 if arithmetic exception occurs, otherwise 0 */
356  )
357 { /*lint --e{715} */
358  SCIP_NLPIPROBLEM* problem;
359  real val;
360  int j;
361 
362  problem = (SCIP_NLPIPROBLEM*)(void*)iuser;
363  assert(problem != NULL);
364 
365  *errflag = 0;
366  for( j = 0; j < *m; ++j )
367  {
368  if( SCIPnlpiOracleEvalConstraintValue(problem->oracle, j, x, &val) != SCIP_OKAY || !SCIPisFinite(val) )
369  {
370  *errflag = 1;
371  SCIPdebugMessage("arithmetic exception in confun for constraint %d\n", j);
372  break;
373  }
374  c[j] = val;
375  }
376 }
377 
378 /** Objective gradient and Jacobian evaluation
379  *
380  * \note If an arithmetic exception occurred, then the gradients must not be modified.
381  */
382 void
383 F77_FUNC(gradient,GRADIENT)(
384  fint* n, /**< number of variables */
385  fint* m, /**< number of constraints */
386  fint* mxa, /**< actual number of entries in a */
387  real* x, /**< value of current variables (array of length n) */
388  real* a, /**< Jacobian matrix entries */
389  fint* la, /**< Jacobian index information: column indices and pointers to start of each row */
390  fint* maxa, /**< maximal size of a */
391  real* user, /**< user real workspace */
392  fint* iuser, /**< user integer workspace */
393  fint* errflag /**< set to 1 if arithmetic exception occurs, otherwise 0 */
394  )
395 { /*lint --e{715} */
396  SCIP_NLPIPROBLEM* problem;
397  SCIP_Real dummy;
398 
399  problem = (SCIP_NLPIPROBLEM*)(void*)iuser;
400  assert(problem != NULL);
401  assert(problem->evalbuffer != NULL);
402  assert(problem->evalbufsize >= *maxa);
403 
404  *errflag = 1;
405 
406  if( SCIPnlpiOracleEvalObjectiveGradient(problem->oracle, x, TRUE, &dummy, problem->evalbuffer) == SCIP_OKAY )
407  {
408  if( SCIPnlpiOracleEvalJacobian(problem->oracle, x, TRUE, NULL, problem->evalbuffer+*n) == SCIP_OKAY )
409  {
410  BMScopyMemoryArray(a, problem->evalbuffer, *maxa);
411  *errflag = 0;
412  }
413  else
414  {
415  SCIPdebugMessage("arithmetic exception in gradient for constraints\n");
416  }
417  }
418  else
419  {
420  SCIPdebugMessage("arithmetic exception in gradient for objective\n");
421  }
422 }
423 
424 /* Objective gradient evaluation */
425 /*
426 void F77_FUNC(objgrad,OBJGRAD)(
427  fint*,
428  fint*,
429  fint*,
430  real*,
431  real*,
432  fint*,
433  fint*,
434  real*,
435  fint*,
436  fint*
437  )
438 */
439 void F77_FUNC(objgrad,OBJGRAD)(void)
440 {
441  SCIPerrorMessage("Objgrad not implemented. Should not be called.\n");
442 }
443 
444 /** Hessian of the Lagrangian evaluation
445  *
446  * phase = 1 : Hessian of the Lagrangian without objective Hessian
447  * phase = 2 : Hessian of the Lagrangian (including objective Hessian)
448  *
449  * \note If an arithmetic exception occurred, then the Hessian must not be modified.
450  */
451 void
452 F77_FUNC(hessian,HESSIAN)(
453  real* x, /**< value of current variables (array of length n) */
454  fint* n, /**< number of variables */
455  fint* m, /**< number of constraints */
456  fint* phase, /**< indicates what kind of Hessian matrix is required */
457  real* lam, /**< Lagrangian multipliers (array of length n+m) */
458  real* ws, /**< real workspace for Hessian, passed to Wdotd */
459  fint* lws, /**< integer workspace for Hessian, passed to Wdotd */
460  real* user, /**< user real workspace */
461  fint* iuser, /**< user integer workspace */
462  fint* l_hess, /**< space of Hessian real storage ws. On entry: maximal space allowed, on exit: actual amount used */
463  fint* li_hess, /**< space of Hessian integer storage lws. On entry: maximal space allowed, on exit: actual amount used */
464  fint* errflag /**< set to 1 if arithmetic exception occurs, otherwise 0 */
465  )
466 { /*lint --e{715} */
467  SCIP_NLPIPROBLEM* problem;
468  SCIP_Real* lambda;
469  int nnz;
470  int i;
471 
472  problem = (SCIP_NLPIPROBLEM*)(void*)iuser;
473  assert(problem != NULL);
474  assert(problem->evalbuffer != NULL);
475 
476  nnz = problem->hessiannz[0]-1;
477  assert(problem->evalbufsize >= nnz);
478 
479  *errflag = 1;
480 
481  /* initialize lambda to -lam */
482  BMSallocMemoryArray(&lambda, *m);
483  for( i = 0; i < *m; ++i )
484  lambda[i] = -lam[*n+i];
485 
486  if( SCIPnlpiOracleEvalHessianLag(problem->oracle, x, TRUE, (*phase == 1) ? 0.0 : 1.0, lambda, problem->evalbuffer) == SCIP_OKAY )
487  {
488  *l_hess = nnz;
489 
490  BMScopyMemoryArray(ws, problem->evalbuffer, nnz);
491 
492  *errflag = 0;
493 
494  /* copy the complete problem->hessiannz into lws */
495  for( i = 0; i < nnz + *n + 2; ++i )
496  lws[i] = problem->hessiannz[i];
497  *li_hess = nnz + *n + 2;
498  }
499  else
500  {
501  SCIPdebugMessage("arithmetic exception in hessian\n");
502  }
503 
504  BMSfreeMemoryArray(&lambda);
505 }
506 
507 
508 
509 static
511  BMS_BLKMEM* blkmem, /**< block memory */
512  SCIP_NLPIORACLE* oracle, /**< NLPI oracle */
513  fint** la, /**< buffer to store pointer to sparsity structure */
514  int* lasize, /**< buffer to store length of *la array */
515  real** a /**< buffer to store pointer to value buffer */
516  )
517 {
518  const int* offset;
519  const int* col;
520  int nnz; /* number of nonzeros in Jacobian */
521  int nvars;
522  int ncons;
523  int i;
524  int c;
525 
526  assert(la != NULL);
527  assert(lasize != NULL);
528  assert(a != NULL);
529  assert(*la == NULL);
530  assert(*a == NULL);
531 
532  nvars = SCIPnlpiOracleGetNVars(oracle);
533  ncons = SCIPnlpiOracleGetNConstraints(oracle);
534 
535  /* get jacobian sparsity in oracle format: offset are rowstarts in col and col are column indices */
536  SCIP_CALL( SCIPnlpiOracleGetJacobianSparsity(oracle, &offset, &col) );
537  nnz = offset[ncons];
538 
539  /* la stores both column indices (first) and rowstarts (at the end) of the objective gradient and Jacobian
540  * For the objective gradient, always n entries are taken, the Jacobian is stored sparse
541  * la(0) = n+nnz+1 position where rowstarts start in la
542  * la(j) column index of objective gradient or Jacobian row, rowwise
543  * la(la(0)) position of first nonzero element for objective gradient in a()
544  * la(la(0)+i) position of first nonzero element for constraint i gradient in a(), i=1..m
545  * la(la(0)+m+1) = n+nnz first unused position in a
546  * where n = nvars and m = ncons
547  */
548 
549  /* need space for la(0) and column indices and rowstarts (1+ncons+1 for objective, constraints, and end (offsets[ncons])) */
550  *lasize = 1 + (nvars+nnz) + 1+ncons + 1;
551  SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, la, *lasize) );
552 
553  (*la)[0] = nvars+nnz+1;
554 
555  /* the objective gradient is stored in sparse form */
556  for( i = 0; i < nvars; ++i )
557  (*la)[1+i] = 1+i; /* shift by 1 for Fortran */
558  (*la)[(*la)[0]] = 1; /* objective entries start at the beginning in a, shift by 1 for Fortran */
559 
560  /* column indicies are as given by col */
561  for( i = 0; i < nnz; ++i )
562  (*la)[1+nvars+i] = col[i] + 1; /* shift by 1 for Fortran */
563 
564  /* rowstarts are as given by offset, plus extra for objective gradient */
565  for( c = 0; c <= ncons; ++c )
566  (*la)[(*la)[0]+1+c] = offset[c] + nvars + 1; /* shift by nvars for objective, shift by 1 for Fortran */
567 
568  SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, a, nvars + nnz) );
569 
570 #if 0 /* enable for debugging Jacobian */
571  for( i = 0; i < 1 + (nvars+nnz) + 1+ncons + 1; ++i )
572  printf("la[%2d] = %2d\n", i, (*la)[i]);
573 #endif
574 
575  return SCIP_OKAY;
576 }
577 
578 static
580  BMS_BLKMEM* blkmem, /**< block memory */
581  SCIP_NLPIORACLE* oracle, /**< NLPI oracle */
582  fint** la, /**< buffer to store pointer to Hessian sparsity structure */
583  int* lasize /**< buffer to store length of *la array */
584  )
585 {
586  const int* offset;
587  const int* col;
588  int nnz; /* number of nonzeros in Jacobian */
589  int nvars;
590  int i;
591  int v;
592 
593  assert(la != NULL);
594  assert(lasize != NULL);
595  assert(*la == NULL);
596 
597  nvars = SCIPnlpiOracleGetNVars(oracle);
598 
599  /* get Hessian sparsity in oracle format: offset are rowstarts in col and col are column indices */
600  SCIP_CALL( SCIPnlpiOracleGetHessianLagSparsity(oracle, &offset, &col) );
601  nnz = offset[nvars];
602 
603  /* la stores both column indices (first) and rowstarts (at the end) of the (sparse) Hessian
604  * la(0) = nnz+1 position where rowstarts start in la
605  * la(j) column index of Hessian row, rowwise
606  * la(la(0)+i) position of first nonzero element for row i, i = 0..n-1
607  * la(la(0)+n) = nnz first unused position in Hessian
608  * where n = nvars
609  */
610 
611  /* need space for la(0) and column indices and rowstarts
612  * 1 for la(0)
613  * nnz for column indices
614  * nvars for rowstarts
615  * 1 for first unused position
616  */
617  *lasize = 1 + nnz + nvars + 1;
618  SCIP_ALLOC( BMSallocBlockMemoryArray(blkmem, la, *lasize) );
619 
620  (*la)[0] = nnz+1;
621 
622  /* column indicies are as given by col */
623  for( i = 0; i < nnz; ++i )
624  (*la)[1+i] = col[i] + 1; /* shift by 1 for Fortran */
625 
626  /* rowstarts are as given by offset */
627  for( v = 0; v <= nvars; ++v )
628  (*la)[(*la)[0]+v] = offset[v] + 1; /* shift by 1 for Fortran */
629 
630  F77_FUNC(hessc,HESSC).phl = 1;
631 
632 #if 0 /* enable for debugging Hessian */
633  for( i = 0; i < 1 + nnz + nvars + 1; ++i )
634  printf("lw[%2d] = %2d\n", i, (*la)[i]);
635 #endif
636 
637  return SCIP_OKAY;
638 }
639 
640 /** setup starting point for FilterSQP */
641 static
643  SCIP_NLPIDATA* data, /**< NLPI data */
644  SCIP_NLPIPROBLEM* problem, /**< NLPI problem */
645  real* x, /**< array to store initial values */
646  SCIP_Bool* success /**< whether we could setup a point in which functions could be evaluated */
647  )
648 {
649  int i;
650  int n;
651  SCIP_Real val;
652 
653  assert(data != NULL);
654  assert(problem != NULL);
655  assert(x != NULL);
656  assert(success != NULL);
657 
658  n = SCIPnlpiOracleGetNVars(problem->oracle);
659 
660  /* setup starting point */
661  if( problem->initguess != NULL )
662  {
663  for( i = 0; i < n; ++i )
664  x[i] = problem->initguess[i];
665  }
666  else
667  {
668  SCIP_Real lb;
669  SCIP_Real ub;
670 
671  SCIPdebugMessage("FilterSQP started without initial primal values; make up something by projecting 0 onto variable bounds and perturb\n");
672 
673  if( data->randnumgen == NULL )
674  {
675  SCIP_CALL( SCIPrandomCreate(&data->randnumgen, data->blkmem, RANDSEED) );
676  }
677 
678  for( i = 0; i < n; ++i )
679  {
680  lb = SCIPnlpiOracleGetVarLbs(problem->oracle)[i];
681  ub = SCIPnlpiOracleGetVarUbs(problem->oracle)[i];
682  if( lb > 0.0 )
683  x[i] = SCIPrandomGetReal(data->randnumgen, lb, lb + MAXPERTURB*MIN(1.0, ub-lb));
684  else if( ub < 0.0 )
685  x[i] = SCIPrandomGetReal(data->randnumgen, ub - MAXPERTURB*MIN(1.0, ub-lb), ub);
686  else
687  x[i] = SCIPrandomGetReal(data->randnumgen, MAX(lb, -MAXPERTURB*MIN(1.0, ub-lb)), MIN(ub, MAXPERTURB*MIN(1.0, ub-lb)));
688  }
689  }
690 
691  /* check whether objective and constraints can be evaluated and differentiated once in starting point
692  * NOTE: this does not check whether hessian can be computed!
693  */
694  *success = SCIPnlpiOracleEvalObjectiveValue(problem->oracle, x, &val) == SCIP_OKAY && SCIPisFinite(val);
695  *success &= SCIPnlpiOracleEvalObjectiveGradient(problem->oracle, x, FALSE, &val, problem->evalbuffer) == SCIP_OKAY; /*lint !e514*/
696  i = 0;
697  for( ; *success && i < SCIPnlpiOracleGetNConstraints(problem->oracle); ++i )
698  *success = SCIPnlpiOracleEvalConstraintValue(problem->oracle, i, x, &val) == SCIP_OKAY && SCIPisFinite(val);
699  *success &= SCIPnlpiOracleEvalJacobian(problem->oracle, x, FALSE, NULL, problem->evalbuffer) == SCIP_OKAY; /*lint !e514*/
700 
701  if( !*success )
702  {
703  SCIPdebugMessage("could not evaluate or constraint %d in %s starting point or Jacobian not available\n", i-1, problem->initguess != NULL ? "provided" : "made up");
704 
705  if( problem->initguess != NULL )
706  {
707  /* forget given starting point and try to make up our own */
708  BMSfreeBlockMemoryArray(data->blkmem, &problem->initguess, problem->varssize);
709  SCIP_CALL( setupStart(data, problem, x, success) );
710  }
711  }
712 
713  return SCIP_OKAY;
714 }
715 
716 
717 
718 /** sets the solstat and termstat to unknown and other, resp. */
719 static
721  SCIP_NLPIPROBLEM* problem /**< data structure of problem */
722  )
723 {
724  assert(problem != NULL);
725 
726  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
727  problem->termstat = SCIP_NLPTERMSTAT_OTHER;
728 }
729 
730 
731 /** processes results from FilterSQP call */
732 static
734  SCIP_NLPIDATA* nlpidata, /**< NLPI data */
735  SCIP_NLPIPROBLEM* problem, /**< NLPI problem */
736  fint ifail, /**< fail flag from FilterSQP call */
737  real* x, /**< primal solution values from FilterSQP call, or NULL if stopped before filtersqp got called */
738  real* lam /**< dual solution values from FilterSQP call, or NULL if stopped before filtersqp got called */
739  )
740 {
741  int i;
742  int nvars;
743  int ncons;
744 
745  assert(problem != NULL);
746  assert(ifail >= 0);
747  assert((x != NULL) == (lam != NULL));
748 
749  problem->solvetime = timeelapsed(nlpidata);
750 
751  nvars = SCIPnlpiOracleGetNVars(problem->oracle);
752  ncons = SCIPnlpiOracleGetNConstraints(problem->oracle);
753 
754  if( ifail <= 8 && x != NULL )
755  {
756  /* FilterSQP terminated somewhat normally -> store away solution */
757 
758  /* make sure we have memory for solution */
759  if( problem->primalvalues == NULL )
760  {
761  assert(problem->varssize >= nvars); /* ensured in nlpiAddVariables */
762  assert(problem->conssize >= ncons); /* ensured in nlpiAddConstraints */
763  assert(problem->consdualvalues == NULL); /* if primalvalues == NULL, then also consdualvalues should be NULL */
764  assert(problem->varlbdualvalues == NULL); /* if primalvalues == NULL, then also varlbdualvalues should be NULL */
765  assert(problem->varubdualvalues == NULL); /* if primalvalues == NULL, then also varubdualvalues should be NULL */
766 
767  SCIP_ALLOC( BMSallocBlockMemoryArray(nlpidata->blkmem, &problem->primalvalues, problem->varssize) );
768  SCIP_ALLOC( BMSallocBlockMemoryArray(nlpidata->blkmem, &problem->consdualvalues, problem->conssize) );
769  SCIP_ALLOC( BMSallocBlockMemoryArray(nlpidata->blkmem, &problem->varlbdualvalues, problem->varssize) );
770  SCIP_ALLOC( BMSallocBlockMemoryArray(nlpidata->blkmem, &problem->varubdualvalues, problem->varssize) );
771  }
772  else
773  {
774  assert(problem->consdualvalues != NULL);
775  assert(problem->varlbdualvalues != NULL);
776  assert(problem->varubdualvalues != NULL);
777  }
778 
779  for( i = 0; i < nvars; ++i )
780  {
781  problem->primalvalues[i] = x[i];
782 
783  /* if dual from FilterSQP is positive, then it belongs to the lower bound, otherwise to the upper bound */
784  problem->varlbdualvalues[i] = MAX(0.0, lam[i]);
785  problem->varubdualvalues[i] = MAX(0.0, -lam[i]);
786  }
787 
788  /* duals from FilterSQP are negated */
789  for( i = 0; i < ncons; ++i )
790  problem->consdualvalues[i] = -lam[nvars + i];
791  }
792 
793  /* translate ifail to solution and termination status and decide whether we could warmstart next */
794  problem->warmstart = FALSE;
795  switch( ifail )
796  {
797  case 0: /* successful run, solution found */
798  assert(problem->rstat[4] <= problem->feastol); /* should be feasible */
799  problem->solstat = (problem->rstat[0] <= problem->opttol ? SCIP_NLPSOLSTAT_LOCOPT : SCIP_NLPSOLSTAT_FEASIBLE);
800  problem->termstat = SCIP_NLPTERMSTAT_OKAY;
801  problem->warmstart = TRUE;
802  break;
803  case 1: /* unbounded, feasible point with f(x) <= fmin */
804  assert(problem->rstat[4] <= problem->feastol); /* should be feasible */
806  if( problem->fmin == DEFAULT_LOBJLIM ) /*lint !e777*/
807  problem->termstat = SCIP_NLPTERMSTAT_OKAY; /* fmin was not set */
808  else
810  break;
811  case 2: /* linear constraints are inconsistent */
813  problem->termstat = SCIP_NLPTERMSTAT_OKAY;
814  break;
815  case 3: /* (locally) nonlinear infeasible, minimal-infeasible solution found */
816  /* problem->solstat = (problem->rstat[0] <= problem->opttol ? SCIP_NLPSOLSTAT_LOCINFEASIBLE : SCIP_NLPSOLSTAT_UNKNOWN); */
817  problem->solstat = SCIP_NLPSOLSTAT_LOCINFEASIBLE; /* TODO FilterSQP does not set rstat[0] in this case, assuming local infeasibility is valid */
818  problem->termstat = SCIP_NLPTERMSTAT_OKAY;
819  problem->warmstart = TRUE;
820  break;
821  case 4: /* terminate at point with h(x) <= eps (constraint violation below epsilon) but QP infeasible */
822  assert(problem->rstat[4] <= problem->feastol); /* should be feasible */
825  problem->warmstart = TRUE;
826  break;
827  case 5: /* termination with rho < eps (trust region radius below epsilon) */
828  if( problem->rstat[4] <= problem->feastol )
830  else
831  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
833  problem->warmstart = TRUE;
834  break;
835  case 6: /* termination with iter > max_iter */
836  if( problem->rstat[4] <= problem->feastol )
838  else
839  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
840  problem->termstat = SCIP_NLPTERMSTAT_ITLIM;
841  problem->warmstart = TRUE;
842  break;
843  case 7: /* crash in user routine (IEEE error) could not be resolved, or timelimit reached */
844  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
845  if( problem->solvetime >= problem->maxtime )
846  {
847  problem->termstat = SCIP_NLPTERMSTAT_TILIM;
848  problem->warmstart = TRUE;
849  }
850  else
852  break;
853  case 8: /* unexpect ifail from QP solver */
854  if( problem->rstat[4] <= problem->feastol )
856  else
857  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
858  problem->termstat = SCIP_NLPTERMSTAT_OTHER;
859  break;
860  case 9: /* not enough REAL workspace */
861  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
863  break;
864  case 10: /* not enough INTEGER workspace */
865  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
867  break;
868  default:
869  problem->solstat = SCIP_NLPSOLSTAT_UNKNOWN;
870  problem->termstat = SCIP_NLPTERMSTAT_OTHER;
871  break;
872  }
873 
874  return SCIP_OKAY;
875 }
876 
877 
878 /** calculate memory size for dynamically allocated arrays (copied from scip/set.c) */
879 static
881  int num /**< minimum number of entries to store */
882  )
883 {
884  int size;
885 
886  /* calculate the size with this loop, such that the resulting numbers are always the same (-> block memory) */
887  size = 4;
888  while( size < num )
889  size = (int)(1.2 * size + 4);
890 
891  return size;
892 }
893 
894 
895 /*
896  * Callback methods of NLP solver interface
897  */
898 
899 /** copy method of NLP interface (called when SCIP copies plugins)
900  *
901  * input:
902  * - blkmem block memory in target SCIP
903  * - sourcenlpi the NLP interface to copy
904  * - targetnlpi buffer to store pointer to copy of NLP interface
905  */
906 static
907 SCIP_DECL_NLPICOPY( nlpiCopyFilterSQP )
908 {
909  SCIP_NLPIDATA* sourcedata;
910 
911  assert(sourcenlpi != NULL);
912  assert(targetnlpi != NULL);
913 
914  sourcedata = SCIPnlpiGetData(sourcenlpi);
915  assert(sourcedata != NULL);
916 
917  SCIP_CALL( SCIPcreateNlpSolverFilterSQP(blkmem, targetnlpi) );
918  assert(*targetnlpi != NULL);
919 
920  SCIP_CALL( SCIPnlpiSetRealPar(*targetnlpi, NULL, SCIP_NLPPAR_INFINITY, sourcedata->infinity) );
921  SCIP_CALL( SCIPnlpiSetMessageHdlr(*targetnlpi, sourcedata->messagehdlr) );
922 
923  return SCIP_OKAY; /*lint !e527*/
924 } /*lint !e715*/
925 
926 /** destructor of NLP interface to free nlpi data
927  *
928  * input:
929  * - nlpi datastructure for solver interface
930  */
931 static
932 SCIP_DECL_NLPIFREE( nlpiFreeFilterSQP )
933 {
934  SCIP_NLPIDATA* data;
935 
936  assert(nlpi != NULL);
937 
938  data = SCIPnlpiGetData(nlpi);
939  assert(data != NULL);
940 
941  if( data->randnumgen != NULL )
942  {
943  SCIPrandomFree(&data->randnumgen, data->blkmem);
944  }
945 
946  BMSfreeBlockMemory(data->blkmem, &data);
947  assert(data == NULL);
948 
949  return SCIP_OKAY;
950 }
951 
952 /** gets pointer for NLP solver
953  *
954  * to do dirty stuff
955  *
956  * input:
957  * - nlpi datastructure for solver interface
958  *
959  * return: void pointer to solver
960  */
961 static
962 SCIP_DECL_NLPIGETSOLVERPOINTER(nlpiGetSolverPointerFilterSQP)
963 {
964  return NULL; /*lint !e527*/
965 } /*lint !e715*/
966 
967 /** creates a problem instance
968  *
969  * input:
970  * - nlpi datastructure for solver interface
971  * - problem pointer to store the problem data
972  * - name name of problem, can be NULL
973  */
974 static
975 SCIP_DECL_NLPICREATEPROBLEM(nlpiCreateProblemFilterSQP)
976 {
977  SCIP_NLPIDATA* data;
978 
979  assert(nlpi != NULL);
980  assert(problem != NULL);
981 
982  data = SCIPnlpiGetData(nlpi);
983  assert(data != NULL);
984 
985  SCIP_ALLOC( BMSallocBlockMemory(data->blkmem, problem) );
986  BMSclearMemory(*problem);
987 
988  SCIP_CALL( SCIPnlpiOracleCreate(data->blkmem, &(*problem)->oracle) );
989  SCIP_CALL( SCIPnlpiOracleSetInfinity((*problem)->oracle, data->infinity) );
990  SCIP_CALL( SCIPnlpiOracleSetProblemName((*problem)->oracle, name) );
991 
992  (*problem)->feastol = DEFAULT_FEASOPTTOL;
993  (*problem)->opttol = DEFAULT_FEASOPTTOL;
994  (*problem)->fmin = DEFAULT_LOBJLIM;
995  (*problem)->maxtime = DBL_MAX;
996  (*problem)->maxiter = DEFAULT_MAXITER;
997  (*problem)->iprint = 0;
998 
999  invalidateSolution(*problem);
1000 
1001  return SCIP_OKAY; /*lint !e527*/
1002 } /*lint !e715*/
1003 
1004 /** free a problem instance
1005  *
1006  * input:
1007  * - nlpi datastructure for solver interface
1008  * - problem pointer where problem data is stored
1009  */
1010 static
1011 SCIP_DECL_NLPIFREEPROBLEM(nlpiFreeProblemFilterSQP)
1012 {
1013  SCIP_NLPIDATA* data;
1014 
1015  assert(nlpi != NULL);
1016  assert(problem != NULL);
1017  assert(*problem != NULL);
1018 
1019  data = SCIPnlpiGetData(nlpi);
1020  assert(data != NULL);
1021 
1022  if( (*problem)->oracle != NULL )
1023  {
1024  SCIP_CALL( SCIPnlpiOracleFree(&(*problem)->oracle) );
1025  }
1026 
1027  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->initguess, (*problem)->varssize);
1028  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->primalvalues, (*problem)->varssize);
1029  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->consdualvalues, (*problem)->conssize);
1030  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->varlbdualvalues, (*problem)->varssize);
1031  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->varubdualvalues, (*problem)->varssize);
1032  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->cstype, (*problem)->conssize);
1033  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->s, (*problem)->varssize + (*problem)->conssize); /*lint !e776 */
1034  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->bu, (*problem)->varssize + (*problem)->conssize); /*lint !e776 */
1035  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->bl, (*problem)->varssize + (*problem)->conssize); /*lint !e776 */
1036  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->x, (*problem)->varssize);
1037  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->c, (*problem)->conssize);
1038  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->lam, (*problem)->varssize + (*problem)->conssize); /*lint !e776 */
1039  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->a, (*problem)->la != NULL ? (*problem)->la[0]-1 : 0);
1040  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->la, (*problem)->lasize);
1041  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->hessiannz, (*problem)->hessiannzsize);
1042  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->lws, (*problem)->mxiwk);
1043  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->ws, (*problem)->mxwk);
1044  BMSfreeBlockMemoryArrayNull(data->blkmem, &(*problem)->evalbuffer, (*problem)->evalbufsize);
1045 
1046  BMSfreeBlockMemory(data->blkmem, problem);
1047  assert(*problem == NULL);
1048 
1049  return SCIP_OKAY;
1050 } /*lint !e715*/
1051 
1052 /** gets pointer to solver-internal problem instance
1053  *
1054  * to do dirty stuff
1055  *
1056  * input:
1057  * - nlpi datastructure for solver interface
1058  * - problem datastructure for problem instance
1059  *
1060  * return: void pointer to problem instance
1061  */
1062 static
1063 SCIP_DECL_NLPIGETPROBLEMPOINTER(nlpiGetProblemPointerFilterSQP)
1064 {
1065  return NULL; /*lint !e527*/
1066 } /*lint !e715*/
1067 
1068 /** add variables
1069  *
1070  * input:
1071  * - nlpi datastructure for solver interface
1072  * - problem datastructure for problem instance
1073  * - nvars number of variables
1074  * - lbs lower bounds of variables, can be NULL if -infinity
1075  * - ubs upper bounds of variables, can be NULL if +infinity
1076  * - varnames names of variables, can be NULL
1077  */
1078 static
1079 SCIP_DECL_NLPIADDVARS( nlpiAddVarsFilterSQP )
1080 {
1081  SCIP_NLPIDATA* nlpidata;
1082  int oldnvars;
1083 
1084  assert(nlpi != NULL);
1085  assert(problem != NULL);
1086  assert(problem->oracle != NULL);
1087 
1088  nlpidata = SCIPnlpiGetData(nlpi);
1089  assert(nlpidata != NULL);
1090 
1091  oldnvars = SCIPnlpiOracleGetNVars(problem->oracle);
1092 
1093  SCIP_CALL( SCIPnlpiOracleAddVars(problem->oracle, nvars, lbs, ubs, varnames) );
1094 
1095  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->initguess, problem->varssize);
1096  invalidateSolution(problem);
1097  problem->warmstart = FALSE;
1098 
1099  /* increase variables-related arrays in problem, if necessary */
1100  if( problem->varssize < SCIPnlpiOracleGetNVars(problem->oracle) )
1101  {
1102  int newsize = calcGrowSize(SCIPnlpiOracleGetNVars(problem->oracle));
1103 
1104  if( problem->primalvalues != NULL )
1105  {
1106  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->primalvalues, problem->varssize, newsize) );
1107  }
1108 
1109  if( problem->varlbdualvalues != NULL )
1110  {
1111  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->varlbdualvalues, problem->varssize, newsize) );
1112  }
1113 
1114  if( problem->varubdualvalues != NULL )
1115  {
1116  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->varubdualvalues, problem->varssize, newsize) );
1117  }
1118 
1119  if( problem->x != NULL )
1120  {
1121  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->x, problem->varssize, newsize) );
1122  }
1123 
1124  if( problem->lam != NULL )
1125  {
1126  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->lam, problem->varssize + problem->conssize, newsize + problem->conssize) ); /*lint !e776*/
1127  }
1128 
1129  if( problem->bl != NULL )
1130  {
1131  assert(problem->bu != NULL);
1132  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->bl, problem->varssize + problem->conssize, newsize + problem->conssize) ); /*lint !e776*/
1133  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->bu, problem->varssize + problem->conssize, newsize + problem->conssize) ); /*lint !e776*/
1134  }
1135 
1136  if( problem->s != NULL )
1137  {
1138  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->s, problem->varssize + problem->conssize, newsize + problem->conssize) ); /*lint !e776*/
1139  }
1140 
1141  problem->varssize = newsize;
1142  }
1143 
1144  /* update variable bounds in FilterSQP data */
1145  if( problem->bl != NULL )
1146  {
1147  int i;
1148  int nconss;
1149 
1150  nconss = SCIPnlpiOracleGetNConstraints(problem->oracle);
1151 
1152  /* bl and bu have first variable bounds and then constraint sides
1153  * copy the constraint sides to their new position before putting in the new variable bounds
1154  */
1155  for( i = nconss-1; i >= 0; --i )
1156  {
1157  problem->bl[oldnvars+nvars+i] = problem->bl[oldnvars+i];
1158  problem->bu[oldnvars+nvars+i] = problem->bu[oldnvars+i];
1159  }
1160 
1161  /* set bounds of new variable */
1162  for( i = 0; i < nvars; ++i )
1163  {
1164  problem->bl[oldnvars+i] = lbs[i];
1165  problem->bu[oldnvars+i] = ubs[i];
1166  }
1167  }
1168 
1169  /* gradients information is out of date now (objective gradient is stored in dense form) */
1170  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1171  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1172 
1173  /* Hessian information is out of date now (no new entries in Hessian, but also empty cols shows up in sparsity info) */
1174  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1175 
1176  return SCIP_OKAY;
1177 } /*lint !e715*/
1178 
1179 
1180 /** add constraints
1181  * quadratic coefficients: row oriented matrix for each constraint
1182  *
1183  * input:
1184  * - nlpi datastructure for solver interface
1185  * - problem datastructure for problem instance
1186  * - ncons number of added constraints
1187  * - lhss left hand sides of constraints
1188  * - rhss right hand sides of constraints
1189  * - nlininds number of linear coefficients for each constraint
1190  * may be NULL in case of no linear part
1191  * - lininds indices of variables for linear coefficients for each constraint
1192  * may be NULL in case of no linear part
1193  * - linvals values of linear coefficient for each constraint
1194  * may be NULL in case of no linear part
1195  * - nquadrows number of columns in matrix of quadratic part for each constraint
1196  * may be NULL in case of no quadratic part in any constraint
1197  * - quadrowidxs indices of variables for which a quadratic part is specified
1198  * may be NULL in case of no quadratic part in any constraint
1199  * - quadoffsets start index of each rows quadratic coefficients in quadinds[.] and quadvals[.]
1200  * indices are given w.r.t. quadrowidxs., i.e., quadoffsets[.][i] gives the start index of row quadrowidxs[.][i] in quadvals[.]
1201  * quadoffsets[.][nquadrows[.]] gives length of quadinds[.] and quadvals[.]
1202  * entry of array may be NULL in case of no quadratic part
1203  * may be NULL in case of no quadratic part in any constraint
1204  * - quadinds column indices w.r.t. quadrowidxs, i.e., quadrowidxs[quadinds[.][i]] gives the index of the variable corresponding
1205  * to entry i, entry of array may be NULL in case of no quadratic part
1206  * may be NULL in case of no quadratic part in any constraint
1207  * - quadvals coefficient values
1208  * entry of array may be NULL in case of no quadratic part
1209  * may be NULL in case of no quadratic part in any constraint
1210  * - exprvaridxs indices of variables in expression tree, maps variable indices in expression tree to indices in nlp
1211  * entry of array may be NULL in case of no expression tree
1212  * may be NULL in case of no expression tree in any constraint
1213  * - exprtrees expression tree for nonquadratic part of constraints
1214  * entry of array may be NULL in case of no nonquadratic part
1215  * may be NULL in case of no nonquadratic part in any constraint
1216  * - names of constraints, may be NULL or entries may be NULL
1217  */
1218 static
1219 SCIP_DECL_NLPIADDCONSTRAINTS( nlpiAddConstraintsFilterSQP )
1220 {
1221  SCIP_NLPIDATA* nlpidata;
1222  int oldnconss;
1223 
1224  assert(nlpi != NULL);
1225  assert(problem != NULL);
1226  assert(problem->oracle != NULL);
1227 
1228  nlpidata = SCIPnlpiGetData(nlpi);
1229  assert(nlpidata != NULL);
1230 
1231  oldnconss = SCIPnlpiOracleGetNConstraints(problem->oracle);
1232 
1233  SCIP_CALL( SCIPnlpiOracleAddConstraints(problem->oracle,
1234  ncons, lhss, rhss,
1235  nlininds, lininds, linvals,
1236  nquadelems, quadelems,
1237  exprvaridxs, exprtrees, names) );
1238 
1239  invalidateSolution(problem);
1240  problem->warmstart = FALSE;
1241 
1242  /* increase constraints-related arrays in problem, if necessary */
1243  if( SCIPnlpiOracleGetNConstraints(problem->oracle) > problem->conssize )
1244  {
1245  int newsize = calcGrowSize(SCIPnlpiOracleGetNConstraints(problem->oracle));
1246 
1247  if( problem->consdualvalues != NULL )
1248  {
1249  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->consdualvalues, problem->conssize, newsize) );
1250  }
1251 
1252  if( problem->c != NULL )
1253  {
1254  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->c, problem->conssize, newsize) );
1255  }
1256 
1257  if( problem->lam != NULL )
1258  {
1259  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->lam, problem->varssize + problem->conssize, problem->varssize + newsize) ); /*lint !e776*/
1260  }
1261 
1262  if( problem->bl != NULL )
1263  {
1264  assert(problem->bu != NULL);
1265  assert(problem->cstype != NULL);
1266  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->bl, problem->varssize + problem->conssize, problem->varssize + newsize) ); /*lint !e776*/
1267  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->bu, problem->varssize + problem->conssize, problem->varssize + newsize) ); /*lint !e776*/
1268  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->cstype, problem->conssize, newsize) );
1269  }
1270 
1271  if( problem->s != NULL )
1272  {
1273  SCIP_ALLOC( BMSreallocBlockMemoryArray(nlpidata->blkmem, &problem->s, problem->varssize + problem->conssize, problem->varssize + newsize) ); /*lint !e776*/
1274  }
1275 
1276  problem->conssize = newsize;
1277  }
1278 
1279  /* update constraint sides and type in FilterSQP data */
1280  if( problem->bl != NULL )
1281  {
1282  int i;
1283  int nvars;
1284 
1285  nvars = SCIPnlpiOracleGetNVars(problem->oracle);
1286 
1287  for( i = 0; i < ncons; ++i )
1288  {
1289  problem->bl[nvars+oldnconss+i] = lhss[i];
1290  problem->bu[nvars+oldnconss+i] = rhss[i];
1291  problem->cstype[oldnconss+i] = SCIPnlpiOracleGetConstraintDegree(problem->oracle, oldnconss+i) <= 1 ? 'L' : 'N';
1292  }
1293  }
1294 
1295  /* gradients information is out of date now */
1296  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1297  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1298 
1299  /* Hessian information is out of date now */
1300  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1301 
1302  return SCIP_OKAY;
1303 } /*lint !e715*/
1304 
1305 /** sets or overwrites objective, a minimization problem is expected
1306  * May change sparsity pattern.
1307  *
1308  * input:
1309  * - nlpi datastructure for solver interface
1310  * - problem datastructure for problem instance
1311  * - nlins number of linear variables
1312  * - lininds variable indices
1313  * may be NULL in case of no linear part
1314  * - linvals coefficient values
1315  * may be NULL in case of no linear part
1316  * - nquadcols number of columns in matrix of quadratic part
1317  * - quadcols indices of variables for which a quadratic part is specified
1318  * may be NULL in case of no quadratic part
1319  * - quadoffsets start index of each rows quadratic coefficients in quadinds and quadvals
1320  * quadoffsets[.][nquadcols] gives length of quadinds and quadvals
1321  * may be NULL in case of no quadratic part
1322  * - quadinds column indices
1323  * may be NULL in case of no quadratic part
1324  * - quadvals coefficient values
1325  * may be NULL in case of no quadratic part
1326  * - exprvaridxs indices of variables in expression tree, maps variable indices in expression tree to indices in nlp
1327  * may be NULL in case of no expression tree
1328  * - exprtree expression tree for nonquadratic part of objective function
1329  * may be NULL in case of no nonquadratic part
1330  * - constant objective value offset
1331  */
1332 static
1333 SCIP_DECL_NLPISETOBJECTIVE( nlpiSetObjectiveFilterSQP )
1334 {
1335  SCIP_NLPIDATA* nlpidata;
1336 
1337  assert(nlpi != NULL);
1338  assert(problem != NULL);
1339  assert(problem->oracle != NULL);
1340 
1341  nlpidata = SCIPnlpiGetData(nlpi);
1342  assert(nlpidata != NULL);
1343 
1344  SCIP_CALL( SCIPnlpiOracleSetObjective(problem->oracle,
1345  constant, nlins, lininds, linvals,
1346  nquadelems, quadelems,
1347  exprvaridxs, exprtree) );
1348 
1349  invalidateSolution(problem);
1350 
1351  /* gradients info (la,a) should still be ok, as objective gradient is stored in dense form */
1352 
1353  /* Hessian information is out of date now */
1354  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1355 
1356  return SCIP_OKAY;
1357 } /*lint !e715*/
1358 
1359 /** change variable bounds
1360  *
1361  * input:
1362  * - nlpi datastructure for solver interface
1363  * - problem datastructure for problem instance
1364  * - nvars number of variables to change bounds
1365  * - indices indices of variables to change bounds
1366  * - lbs new lower bounds
1367  * - ubs new upper bounds
1368  */
1369 static
1370 SCIP_DECL_NLPICHGVARBOUNDS( nlpiChgVarBoundsFilterSQP )
1371 {
1372  assert(nlpi != NULL);
1373  assert(problem != NULL);
1374  assert(problem->oracle != NULL);
1375 
1376  SCIP_CALL( SCIPnlpiOracleChgVarBounds(problem->oracle, nvars, indices, lbs, ubs) );
1377 
1378  invalidateSolution(problem);
1379 
1380  /* update bounds in FilterSQP data */
1381  if( problem->bl != NULL )
1382  {
1383  int i;
1384 
1385  for( i = 0; i < nvars; ++i )
1386  {
1387  problem->bl[indices[i]] = lbs[i];
1388  problem->bu[indices[i]] = ubs[i];
1389  }
1390  }
1391 
1392  return SCIP_OKAY;
1393 } /*lint !e715*/
1394 
1395 /** change constraint bounds
1396  *
1397  * input:
1398  * - nlpi datastructure for solver interface
1399  * - problem datastructure for problem instance
1400  * - nconss number of constraints to change sides
1401  * - indices indices of constraints to change sides
1402  * - lhss new left hand sides
1403  * - rhss new right hand sides
1404  */
1405 static
1406 SCIP_DECL_NLPICHGCONSSIDES( nlpiChgConsSidesFilterSQP )
1407 {
1408  assert(nlpi != NULL);
1409  assert(problem != NULL);
1410  assert(problem->oracle != NULL);
1411 
1412  SCIP_CALL( SCIPnlpiOracleChgConsSides(problem->oracle, nconss, indices, lhss, rhss) );
1413 
1414  invalidateSolution(problem);
1415 
1416  /* update constraint sense in FilterSQP data */
1417  if( problem->bl != NULL )
1418  {
1419  int i;
1420  int nvars;
1421 
1422  nvars = SCIPnlpiOracleGetNVars(problem->oracle);
1423 
1424  for( i = 0; i < nconss; ++i )
1425  {
1426  problem->bl[nvars+indices[i]] = lhss[i];
1427  problem->bu[nvars+indices[i]] = rhss[i];
1428  }
1429  }
1430 
1431  return SCIP_OKAY;
1432 } /*lint !e715*/
1433 
1434 /** delete a set of variables
1435  *
1436  * input:
1437  * - nlpi datastructure for solver interface
1438  * - problem datastructure for problem instance
1439  * - dstats deletion status of vars; 1 if var should be deleted, 0 if not
1440  *
1441  * output:
1442  * - dstats new position of var, -1 if var was deleted
1443  */
1444 static
1445 SCIP_DECL_NLPIDELVARSET( nlpiDelVarSetFilterSQP )
1446 {
1447  SCIP_NLPIDATA* nlpidata;
1448 
1449  assert(nlpi != NULL);
1450  assert(problem != NULL);
1451  assert(problem->oracle != NULL);
1452 
1453  nlpidata = SCIPnlpiGetData(nlpi);
1454  assert(nlpidata != NULL);
1455 
1456  SCIP_CALL( SCIPnlpiOracleDelVarSet(problem->oracle, dstats) );
1457 
1458  /* @TODO keep initguess and bl, bu for remaining variables? */
1459 
1460  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->initguess, problem->varssize);
1461  invalidateSolution(problem);
1462  problem->warmstart = FALSE;
1463 
1464  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->bl, problem->varssize + problem->conssize); /*lint !e776 */
1465  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->bu, problem->varssize + problem->conssize); /*lint !e776 */
1466  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->cstype, problem->conssize); /* because we assume that cstype is allocated iff bl is allocated */
1467 
1468  /* gradients information is out of date now (objective gradient is stored in dense form) */
1469  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1470  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1471 
1472  /* Hessian information is out of date now */
1473  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1474 
1475  return SCIP_OKAY;
1476 } /*lint !e715*/
1477 
1478 /** delete a set of constraints
1479  *
1480  * input:
1481  * - nlpi datastructure for solver interface
1482  * - problem datastructure for problem instance
1483  * - dstats deletion status of rows; 1 if row should be deleted, 0 if not
1484  *
1485  * output:
1486  * - dstats new position of row, -1 if row was deleted
1487  */
1488 static
1489 SCIP_DECL_NLPIDELCONSSET( nlpiDelConstraintSetFilterSQP )
1490 {
1491  SCIP_NLPIDATA* nlpidata;
1492 
1493  assert(nlpi != NULL);
1494  assert(problem != NULL);
1495  assert(problem->oracle != NULL);
1496 
1497  nlpidata = SCIPnlpiGetData(nlpi);
1498  assert(nlpidata != NULL);
1499 
1500  SCIP_CALL( SCIPnlpiOracleDelConsSet(problem->oracle, dstats) );
1501 
1502  invalidateSolution(problem);
1503  problem->warmstart = FALSE;
1504 
1505  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->bl, problem->varssize + problem->conssize); /*lint !e776 */
1506  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->bu, problem->varssize + problem->conssize); /*lint !e776 */
1507  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->cstype, problem->conssize); /* because we assume that cstype is allocated iff bl is allocated */
1508 
1509  /* gradients information is out of date now */
1510  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1511  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1512 
1513  /* Hessian information is out of date now */
1514  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1515 
1516  return SCIP_OKAY;
1517 } /*lint !e715*/
1518 
1519 /** changes (or adds) linear coefficients in a constraint or objective
1520  *
1521  * input:
1522  * - nlpi datastructure for solver interface
1523  * - problem datastructure for problem instance
1524  * - idx index of constraint or -1 for objective
1525  * - nvals number of values in linear constraint to change
1526  * - varidxs indices of variables which coefficient to change
1527  * - vals new values for coefficients
1528  */
1529 static
1530 SCIP_DECL_NLPICHGLINEARCOEFS( nlpiChgLinearCoefsFilterSQP )
1531 {
1532  SCIP_NLPIDATA* nlpidata;
1533 
1534  assert(nlpi != NULL);
1535  assert(problem != NULL);
1536  assert(problem->oracle != NULL);
1537 
1538  nlpidata = SCIPnlpiGetData(nlpi);
1539  assert(nlpidata != NULL);
1540 
1541  SCIP_CALL( SCIPnlpiOracleChgLinearCoefs(problem->oracle, idx, nvals, varidxs, vals) );
1542 
1543  invalidateSolution(problem);
1544 
1545  /* gradients information (la,a) may have changed if elements were added or removed
1546  * (we only care that sparsity doesn't change, not about actual values in a)
1547  * TODO free only if coefficients were added or removed (SCIPnlpiOracleChgLinearCoefs() could give feedback)
1548  */
1549  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1550  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1551 
1552  /* Hessian information should still be ok */
1553 
1554  return SCIP_OKAY;
1555 } /*lint !e715*/
1556 
1557 /** changes (or adds) coefficients in the quadratic part of a constraint or objective
1558  *
1559  * input:
1560  * - nlpi datastructure for solver interface
1561  * - problem datastructure for problem instance
1562  * - idx index of constraint or -1 for objective
1563  * - nentries number of entries in quadratic matrix to change
1564  * - rows row indices of entries in quadratic matrix where values should be changed
1565  * - cols column indices of entries in quadratic matrix where values should be changed
1566  * - values new values for entries in quadratic matrix
1567  */
1568 static
1569 SCIP_DECL_NLPICHGQUADCOEFS( nlpiChgQuadraticCoefsFilterSQP )
1570 {
1571  SCIP_NLPIDATA* nlpidata;
1572 
1573  assert(nlpi != NULL);
1574  assert(problem != NULL);
1575  assert(problem->oracle != NULL);
1576 
1577  nlpidata = SCIPnlpiGetData(nlpi);
1578  assert(nlpidata != NULL);
1579 
1580  SCIP_CALL( SCIPnlpiOracleChgQuadCoefs(problem->oracle, idx, nquadelems, quadelems) );
1581 
1582  invalidateSolution(problem);
1583 
1584  /* update constraint linearity in FilterSQP data, as we might have changed from linear to nonlinear now */
1585  if( problem->cstype != NULL && idx >= 0 )
1586  problem->cstype[idx] = (SCIPnlpiOracleGetConstraintDegree(problem->oracle, idx) <= 1 ? 'L' : 'N');
1587 
1588  /* gradients information (la,a) may have changed if elements were added or removed
1589  * (we only care that sparsity doesn't change, not about actual values in a)
1590  * TODO free only if coefficients were added or removed (SCIPnlpiOracleChgLinearCoefs() could give feedback)
1591  */
1592  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1593  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1594 
1595  /* Hessian sparsity may have changed if elements were added or removed
1596  * TODO free only if coefficients were added or removed (SCIPnlpiOracleChgLinearCoefs() could give feedback)
1597  */
1598  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1599 
1600  return SCIP_OKAY;
1601 } /*lint !e715*/
1602 
1603 /** replaces the expression tree of a constraint or objective
1604  *
1605  * input:
1606  * - nlpi datastructure for solver interface
1607  * - problem datastructure for problem instance
1608  * - idxcons index of constraint or -1 for objective
1609  * - exprvaridxs indices of variables in expression tree, maps variable indices in expression tree to indices in nlp, or NULL
1610  * - exprtree new expression tree for constraint or objective, or NULL to only remove previous tree
1611  */
1612 static
1613 SCIP_DECL_NLPICHGEXPRTREE( nlpiChgExprtreeFilterSQP )
1614 {
1615  SCIP_NLPIDATA* nlpidata;
1616 
1617  assert(nlpi != NULL);
1618  assert(problem != NULL);
1619  assert(problem->oracle != NULL);
1620 
1621  nlpidata = SCIPnlpiGetData(nlpi);
1622  assert(nlpidata != NULL);
1623 
1624  SCIP_CALL( SCIPnlpiOracleChgExprtree(problem->oracle, idxcons, exprvaridxs, exprtree) );
1625 
1626  invalidateSolution(problem);
1627 
1628  /* update constraint linearity in FilterSQP data, as we might have changed from linear to nonlinear now */
1629  if( problem->cstype != NULL && idxcons >= 0 )
1630  problem->cstype[idxcons] = (SCIPnlpiOracleGetConstraintDegree(problem->oracle, idxcons) <= 1 ? 'L' : 'N');
1631 
1632  /* gradients information (la,a) may have changed */
1633  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1634  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1635 
1636  /* Hessian information may have changed */
1637  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1638 
1639  return SCIP_OKAY;
1640 } /*lint !e715*/
1641 
1642 /** change one coefficient in the nonlinear part
1643  *
1644  * input:
1645  * - nlpi datastructure for solver interface
1646  * - problem datastructure for problem instance
1647  * - idxcons index of constraint or -1 for objective
1648  * - idxparam index of parameter
1649  * - value new value for nonlinear parameter
1650  *
1651  * return: Error if parameter does not exist
1652  */
1653 static
1654 SCIP_DECL_NLPICHGNONLINCOEF( nlpiChgNonlinCoefFilterSQP )
1655 {
1656  SCIP_NLPIDATA* nlpidata;
1657 
1658  assert(nlpi != NULL);
1659  assert(problem != NULL);
1660  assert(problem->oracle != NULL);
1661 
1662  nlpidata = SCIPnlpiGetData(nlpi);
1663  assert(nlpidata != NULL);
1664 
1665  SCIP_CALL( SCIPnlpiOracleChgExprParam(problem->oracle, idxcons, idxparam, value) );
1666 
1667  invalidateSolution(problem);
1668 
1669  /* gradients information (la,a) may have changed (?) */
1670  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->a, problem->la != NULL ? problem->la[0]-1 : 0);
1671  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->la, problem->lasize);
1672 
1673  /* Hessian information may have changed (?) */
1674  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->hessiannz, problem->hessiannzsize);
1675 
1676  return SCIP_OKAY;
1677 } /*lint !e715*/
1678 
1679 /** change the constant offset in the objective
1680  *
1681  * input:
1682  * - nlpi datastructure for solver interface
1683  * - problem datastructure for problem instance
1684  * - objconstant new value for objective constant
1685  */
1686 static
1687 SCIP_DECL_NLPICHGOBJCONSTANT( nlpiChgObjConstantFilterSQP )
1688 {
1689  assert(nlpi != NULL);
1690  assert(problem != NULL);
1691  assert(problem->oracle != NULL);
1692 
1693  SCIP_CALL( SCIPnlpiOracleChgObjConstant(problem->oracle, objconstant) );
1694 
1695  invalidateSolution(problem);
1696 
1697  return SCIP_OKAY;
1698 } /*lint !e715*/
1699 
1700 /** sets initial guess for primal variables
1701  *
1702  * input:
1703  * - nlpi datastructure for solver interface
1704  * - problem datastructure for problem instance
1705  * - primalvalues initial primal values for variables, or NULL to clear previous values
1706  * - consdualvalues initial dual values for constraints, or NULL to clear previous values
1707  * - varlbdualvalues initial dual values for variable lower bounds, or NULL to clear previous values
1708  * - varubdualvalues initial dual values for variable upper bounds, or NULL to clear previous values
1709  */
1710 static
1711 SCIP_DECL_NLPISETINITIALGUESS( nlpiSetInitialGuessFilterSQP )
1712 {
1713  SCIP_NLPIDATA* nlpidata;
1714 
1715  assert(nlpi != NULL);
1716  assert(problem != NULL);
1717  assert(problem->oracle != NULL);
1718 
1719  nlpidata = SCIPnlpiGetData(nlpi);
1720  assert(nlpidata != NULL);
1721 
1722  if( primalvalues != NULL )
1723  {
1724  if( problem->initguess == NULL )
1725  {
1726  SCIP_ALLOC( BMSallocBlockMemoryArray(nlpidata->blkmem, &problem->initguess, problem->varssize) );
1727  }
1728  assert(SCIPnlpiOracleGetNVars(problem->oracle) <= problem->varssize);
1729  BMScopyMemoryArray(problem->initguess, primalvalues, SCIPnlpiOracleGetNVars(problem->oracle));
1730  }
1731  else
1732  {
1733  BMSfreeBlockMemoryArrayNull(nlpidata->blkmem, &problem->initguess, problem->varssize);
1734  }
1735 
1736  return SCIP_OKAY;
1737 } /*lint !e715*/
1738 
1739 /** tries to solve NLP
1740  *
1741  * input:
1742  * - nlpi datastructure for solver interface
1743  * - problem datastructure for problem instance
1744  */
1745 static
1746 SCIP_DECL_NLPISOLVE( nlpiSolveFilterSQP )
1747 {
1748  SCIP_NLPIDATA* data;
1749  SCIP_Bool success;
1750  fint n;
1751  fint m;
1752  fint kmax;
1753  fint maxa;
1754  fint maxf;
1755  fint mlp;
1756  fint lh1;
1757  fint nout;
1758  fint ifail;
1759  fint maxiter;
1760  real rho;
1761  real f;
1762  real* user;
1763  fint* iuser;
1764  ftnlen cstype_len = 1;
1765  fint minmxwk;
1766  fint minmxiwk;
1767  int nruns;
1768  int i;
1769 
1770  data = SCIPnlpiGetData(nlpi);
1771  assert(data != NULL);
1772 
1773  /* start measuring time */
1774  data->starttime = gettime();
1775 
1776  /* if fromscratch parameter is set, then we will not warmstart */
1777  if( problem->fromscratch )
1778  problem->warmstart = FALSE;
1779 
1780  n = SCIPnlpiOracleGetNVars(problem->oracle);
1781  m = SCIPnlpiOracleGetNConstraints(problem->oracle);
1782  kmax = n; /* maximal nullspace dimension */
1783  maxf = 100; /* maximal filter length */
1784  mlp = 100; /* maximum level of degeneracy */
1785 
1786  /* TODO eventually, the output should be redirected to the message handler,
1787  * but even to just redirect to some other file, we would have to open the output-unit in Fortran
1788  */
1789  nout = 6; /* output to stdout for now */
1790  ifail = problem->warmstart ? -1 : 0; /* -1 for warmstart, otherwise 0 */
1791  rho = 10.0; /* initial trust-region radius */
1792 
1793  user = (real*)data;
1794  iuser = (fint*)problem;
1795  if( problem->warmstart ) /* if warmstart, then need to keep istat[0] */
1796  memset(problem->istat+1, 0, sizeof(problem->istat)-sizeof(*problem->istat));
1797  else
1798  memset(problem->istat, 0, sizeof(problem->istat));
1799  memset(problem->rstat, 0, sizeof(problem->rstat));
1800  problem->niterations = 0;
1801 
1802  if( problem->x == NULL )
1803  {
1804  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->x, problem->varssize) );
1805  }
1806  if( problem->c == NULL )
1807  {
1808  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->c, problem->conssize) );
1809  }
1810  if( problem->lam == NULL )
1811  {
1812  SCIP_ALLOC( BMSallocClearBlockMemoryArray(data->blkmem, &problem->lam, problem->varssize + problem->conssize) ); /*lint !e776 */
1813  }
1814  else
1815  {
1816  BMSclearMemoryArray(problem->lam, problem->varssize + problem->conssize);
1817  }
1818  if( problem->s == NULL )
1819  {
1820  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->s, problem->varssize + problem->conssize) );
1821  }
1822 
1823  if( problem->la == NULL )
1824  {
1825  /* allocate la, a and initialize la for Objective Gradient and Jacobian */
1826  SCIP_CALL( setupGradients(data->blkmem, problem->oracle, &problem->la, &problem->lasize, &problem->a) );
1827  }
1828  /* maximal number entries in a = nvars+nnz */
1829  maxa = problem->la[0]-1;
1830 
1831  if( problem->hessiannz == NULL )
1832  {
1833  /* allocate and initialize problem->hessiannz for Hessian */
1834  SCIP_CALL( setupHessian(data->blkmem, problem->oracle, &problem->hessiannz, &problem->hessiannzsize) );
1835  }
1836 
1837  /* setup variable bounds, constraint sides, and constraint types */
1838  if( problem->bl == NULL )
1839  {
1840  assert(problem->bu == NULL);
1841  assert(problem->cstype == NULL);
1842 
1843  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->bl, problem->varssize + problem->conssize) );
1844  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->bu, problem->varssize + problem->conssize) );
1845  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->cstype, problem->conssize) );
1846 
1847  BMScopyMemoryArray(problem->bl, SCIPnlpiOracleGetVarLbs(problem->oracle), n);
1848  BMScopyMemoryArray(problem->bu, SCIPnlpiOracleGetVarUbs(problem->oracle), n);
1849  for( i = 0; i < m; ++i )
1850  {
1851  problem->bl[n+i] = SCIPnlpiOracleGetConstraintLhs(problem->oracle, i);
1852  problem->bu[n+i] = SCIPnlpiOracleGetConstraintRhs(problem->oracle, i);
1853  problem->cstype[i] = SCIPnlpiOracleGetConstraintDegree(problem->oracle, i) <= 1 ? 'L' : 'N';
1854  }
1855  }
1856 
1857  /* buffer for evaluation results (used in setupStart already) */
1858  if( problem->evalbufsize < MAX3(n, problem->hessiannz[0], maxa) )
1859  {
1860  int newsize = calcGrowSize(MAX3(n, problem->hessiannz[0], maxa));
1861  SCIP_ALLOC( BMSreallocBlockMemoryArray(data->blkmem, &problem->evalbuffer, problem->evalbufsize, newsize) );
1862  problem->evalbufsize = newsize;
1863  }
1864 
1865  /* setup starting point */
1866  SCIP_CALL( setupStart(data, problem, problem->x, &success) );
1867  if( !success )
1868  {
1869  /* FilterSQP would crash if starting point cannot be evaluated, so give up */
1870  SCIP_CALL( processSolveOutcome(data, problem, 7, NULL, NULL) );
1871  return SCIP_OKAY;
1872  }
1873 
1874  /* setup workspace */
1875  /* initial guess of real workspace size */
1876  /* FilterSQP manual: mxwk = 21*n + 8*m + mlp + 8*maxf + kmax*(kmax+9)/2 + nprof, with nprof = 20*n as a good guess */
1877  /* Bonmin: mxwk = 21*n + 8*m + mlp + 8*maxf + lh1 + kmax*(kmax+9)/2 + mxwk0,
1878  * with lh1 = nnz_h+8+2*n+m and mxwk0 = 2000000 (parameter) */
1879  lh1 = problem->hessiannz[0]-1 + 8 + 2*n + m;
1880  minmxwk = 21*n + 8*m + mlp + 8*maxf + lh1 + kmax*(kmax+9)/2 + MAX(20*n,2000);
1881  if( problem->ws == NULL )
1882  {
1883  problem->mxwk = minmxwk;
1884  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->ws, problem->mxwk) );
1885  }
1886  else if( problem->mxwk < minmxwk )
1887  {
1888  int newsize = calcGrowSize(minmxwk);
1889  SCIP_ALLOC( BMSreallocBlockMemoryArray(data->blkmem, &problem->ws, problem->mxwk, newsize) );
1890  problem->mxwk = newsize;
1891  }
1892 
1893  /* initial guess of integer workspace size */
1894  /* FilterSQP manual: mxiwk = 13*n + 4*m + mlp + 100 + kmax */
1895  /* Bonmin: mxiwk = 13*n + 4*m + mlp + lh1 + kmax + 113 + mxiwk0, with mxiwk0 = 500000 (parameter) */
1896  minmxiwk = 13*n + 4*m + mlp + lh1 + 100 + kmax + 113 + MAX(5*n,5000);
1897  if( problem->lws == NULL )
1898  {
1899  assert(!problem->warmstart);
1900 
1901  problem->mxiwk = minmxiwk;
1902  SCIP_ALLOC( BMSallocBlockMemoryArray(data->blkmem, &problem->lws, problem->mxiwk) );
1903  }
1904  else if( problem->mxiwk < minmxiwk && !problem->warmstart ) /* if warmstart, then lws should remain untouched (n and m didn't change anyway) */
1905  {
1906  int newsize = calcGrowSize(minmxiwk);
1907  SCIP_ALLOC( BMSreallocBlockMemoryArray(data->blkmem, &problem->lws, problem->mxiwk, newsize) );
1908  problem->mxiwk = newsize;
1909  }
1910  /* in case of some evalerrors, not clearing ws could lead to valgrind warnings about use of uninitialized memory */
1911  memset(problem->ws, 0, problem->mxwk * sizeof(real));
1912 
1913  /* from here on we are not thread-safe: if intended for multithread use, then protect filtersqp call with mutex
1914  * NOTE: we need to make sure that we do not return from nlpiSolve before unlocking the mutex
1915  */
1916 #ifndef NPARASCIP
1917  pthread_mutex_lock(&filtersqpmutex);
1918 #endif
1919 
1920  /* initialize global variables from filtersqp */
1921  /* FilterSQP eps is tolerance for both feasibility and optimality, and also for trust-region radius, etc. */
1922  F77_FUNC(nlp_eps_inf,NLP_EPS_INF).eps = MIN(problem->feastol, problem->opttol * OPTTOLFACTOR);
1923  F77_FUNC(nlp_eps_inf,NLP_EPS_INF).infty = SCIPnlpiOracleGetInfinity(problem->oracle);
1924  F77_FUNC(ubdc,UBDC).ubd = 100.0;
1925  F77_FUNC(ubdc,UBDC).tt = 1.25;
1926  F77_FUNC(scalec,SCALEC).scale_mode = 0;
1927 
1928  for( nruns = 1; ; ++nruns )
1929  {
1930  maxiter = problem->maxiter - problem->niterations;
1931 
1932  F77_FUNC(filtersqp,FILTERSQP)(
1933  &n, &m, &kmax, &maxa,
1934  &maxf, &mlp, &problem->mxwk, &problem->mxiwk,
1935  &problem->iprint, &nout, &ifail, &rho,
1936  problem->x, problem->c, &f, &problem->fmin, problem->bl,
1937  problem->bu, problem->s, problem->a, problem->la, problem->ws,
1938  problem->lws, problem->lam, problem->cstype, user,
1939  iuser, &maxiter, problem->istat,
1940  problem->rstat, cstype_len);
1941 
1942  problem->niterations += problem->istat[1];
1943 
1944  assert(ifail <= 10);
1945  /* if ifail >= 8 (probably the workspace was too small), then retry with larger workspace
1946  * if ifail == 0 (local optimal), but absolute violation of KKT too large, then retry with small eps
1947  */
1948  if( ifail < 8 && (ifail != 0 || problem->rstat[0] <= problem->opttol) )
1949  break;
1950 
1951  if( problem->iprint > 0 )
1952  {
1953  SCIPmessagePrintInfo(data->messagehdlr, "FilterSQP terminated with status %d in run %d, absolute KKT violation is %g\n", ifail, nruns, problem->rstat[0]);
1954  }
1955 
1956  /* if iteration or time limit exceeded, then don't retry */
1957  if( problem->niterations >= problem->maxiter || timelimitreached(data, problem) )
1958  {
1959  if( problem->iprint > 0 )
1960  {
1961  SCIPmessagePrintInfo(data->messagehdlr, "Time or iteration limit reached, not retrying\n");
1962  }
1963  break;
1964  }
1965 
1966  /* if maximal number of runs reached, then stop */
1967  if( nruns >= MAXNRUNS )
1968  {
1969  if( problem->iprint > 0 )
1970  {
1971  SCIPmessagePrintInfo(data->messagehdlr, "Run limit reached, not retrying\n");
1972  }
1973  break;
1974  }
1975 
1976  if( ifail == 0 )
1977  {
1978  SCIP_Real epsfactor;
1979 
1980  if( F77_FUNC(nlp_eps_inf,NLP_EPS_INF).eps <= MINEPS )
1981  {
1982  if( problem->iprint > 0 )
1983  {
1984  SCIPmessagePrintInfo(data->messagehdlr, "Already reached minimal epsilon, not retrying\n");
1985  }
1986  break;
1987  }
1988 
1989  epsfactor = problem->opttol / problem->rstat[0];
1990  assert(epsfactor < 1.0); /* because of the if's above */
1991  epsfactor *= OPTTOLFACTOR;
1992 
1993  F77_FUNC(nlp_eps_inf,NLP_EPS_INF).eps = MAX(MINEPS, F77_FUNC(nlp_eps_inf,NLP_EPS_INF).eps * epsfactor);
1994  if( problem->iprint > 0 )
1995  {
1996  SCIPmessagePrintInfo(data->messagehdlr, "Continue with eps = %g\n", F77_FUNC(nlp_eps_inf,NLP_EPS_INF).eps);
1997  }
1998  ifail = -1; /* do warmstart */
1999 
2000  continue;
2001  }
2002 
2003  /* increase real workspace, if ifail = 9 (real workspace too small) or ifail = 8 (unexpected ifail from QP solver, often also when workspace too small) */
2004  if( ifail == 8 || ifail == 9 )
2005  {
2006  int newsize = calcGrowSize(WORKSPACEGROWTHFACTOR*problem->mxwk);
2007  if( BMSreallocBlockMemoryArray(data->blkmem, &problem->ws, problem->mxwk, newsize) == NULL )
2008  {
2009  /* realloc failed: give up NLP solve */
2010  problem->mxwk = 0;
2011  break;
2012  }
2013  problem->mxwk = newsize;
2014  }
2015 
2016  /* increase integer workspace, if ifail = 10 (integer workspace too small) or ifail = 8 (unexpected ifail from QP solver, often also when workspace too small) */
2017  if( ifail == 8 || ifail == 10 )
2018  {
2019  int newsize = calcGrowSize(WORKSPACEGROWTHFACTOR*problem->mxiwk);
2020  if( BMSreallocBlockMemoryArray(data->blkmem, &problem->lws, problem->mxiwk, newsize) == NULL )
2021  {
2022  /* realloc failed: give up NLP solve */
2023  problem->mxiwk = 0;
2024  break;
2025  }
2026  problem->mxiwk = newsize;
2027 
2028  /* better don't try warmstart for the next trial; warmstart requires that lws is untouched, does extending count as touching? */
2029  ifail = 0;
2030  }
2031 
2032  /* reset starting point, in case it was overwritten by failed solve (return can only be SCIP_OKAY, because randnumgen must exist already)
2033  * NOTE: If initguess is NULL (no user-given starting point), then this will result in a slightly different starting point as in the previous setupStart() call (random numbers)
2034  * However, as no point was given, it shouldn't matter which point we actually start from.
2035  */
2036  (void) setupStart(data, problem, problem->x, &success);
2037  assert(success);
2038  }
2039 
2040 #ifndef NPARASCIP
2041  pthread_mutex_unlock(&filtersqpmutex);
2042 #endif
2043 
2044  SCIP_CALL( processSolveOutcome(data, problem, ifail, problem->x, problem->lam) );
2045 
2046  return SCIP_OKAY; /*lint !e527*/
2047 } /*lint !e715*/
2048 
2049 /** gives solution status
2050  *
2051  * input:
2052  * - nlpi datastructure for solver interface
2053  * - problem datastructure for problem instance
2054  *
2055  * return: Solution Status
2056  */
2057 static
2058 SCIP_DECL_NLPIGETSOLSTAT( nlpiGetSolstatFilterSQP )
2059 {
2060  assert(problem != NULL);
2061 
2062  return problem->solstat;
2063 } /*lint !e715*/
2064 
2065 /** gives termination reason
2066  *
2067  * input:
2068  * - nlpi datastructure for solver interface
2069  * - problem datastructure for problem instance
2070  *
2071  * return: Termination Status
2072  */
2073 static
2074 SCIP_DECL_NLPIGETTERMSTAT( nlpiGetTermstatFilterSQP )
2075 {
2076  assert(problem != NULL);
2077 
2078  return problem->termstat;
2079 } /*lint !e715*/
2080 
2081 /** gives primal and dual solution values
2082  *
2083  * solver can return NULL in dual values if not available
2084  * but if solver provides dual values for one side of variable bounds, then it must also provide those for the other side
2085  *
2086  * for a ranged constraint, the dual variable is positive if the right hand side is active and negative if the left hand side is active
2087  *
2088  * input:
2089  * - nlpi datastructure for solver interface
2090  * - problem datastructure for problem instance
2091  * - primalvalues buffer to store pointer to array to primal values, or NULL if not needed
2092  * - consdualvalues buffer to store pointer to array to dual values of constraints, or NULL if not needed
2093  * - varlbdualvalues buffer to store pointer to array to dual values of variable lower bounds, or NULL if not needed
2094  * - varubdualvalues buffer to store pointer to array to dual values of variable lower bounds, or NULL if not needed
2095  * - objval buffer store the objective value, or NULL if not needed
2096  */
2097 static
2098 SCIP_DECL_NLPIGETSOLUTION( nlpiGetSolutionFilterSQP )
2099 {
2100  assert(problem != NULL);
2101 
2102  if( primalvalues != NULL )
2103  {
2104  assert(problem->primalvalues != NULL);
2105 
2106  *primalvalues = problem->primalvalues;
2107  }
2108 
2109  if( consdualvalues != NULL )
2110  {
2111  assert(problem->consdualvalues != NULL);
2112 
2113  *consdualvalues = problem->consdualvalues;
2114  }
2115 
2116  if( varlbdualvalues != NULL )
2117  {
2118  assert(problem->varlbdualvalues != NULL);
2119 
2120  *varlbdualvalues = problem->varlbdualvalues;
2121  }
2122 
2123  if( varubdualvalues != NULL )
2124  {
2125  assert(problem->varubdualvalues != NULL);
2126 
2127  *varubdualvalues = problem->varubdualvalues;
2128  }
2129 
2130  if( objval != NULL )
2131  {
2132  if( problem->primalvalues != NULL )
2133  {
2134  /* TODO store last solution value instead of reevaluating the objective function */
2135  SCIP_CALL( SCIPnlpiOracleEvalObjectiveValue(problem->oracle, problem->primalvalues, objval) );
2136  }
2137  else
2138  *objval = SCIP_INVALID;
2139  }
2140 
2141  return SCIP_OKAY; /*lint !e527*/
2142 } /*lint !e715*/
2143 
2144 /** gives solve statistics
2145  *
2146  * input:
2147  * - nlpi datastructure for solver interface
2148  * - problem datastructure for problem instance
2149  * - statistics pointer to store statistics
2150  *
2151  * output:
2152  * - statistics solve statistics
2153  */
2154 static
2155 SCIP_DECL_NLPIGETSTATISTICS( nlpiGetStatisticsFilterSQP )
2156 {
2157  assert(problem != NULL);
2158 
2159  SCIPnlpStatisticsSetNIterations(statistics, problem->niterations);
2160  SCIPnlpStatisticsSetTotalTime(statistics, problem->solvetime);
2161 
2162  return SCIP_OKAY; /*lint !e527*/
2163 } /*lint !e715*/
2164 
2165 /** gives required size of a buffer to store a warmstart object
2166  *
2167  * input:
2168  * - nlpi datastructure for solver interface
2169  * - problem datastructure for problem instance
2170  * - size pointer to store required size for warmstart buffer
2171  *
2172  * output:
2173  * - size required size for warmstart buffer
2174  */
2175 static
2176 SCIP_DECL_NLPIGETWARMSTARTSIZE( nlpiGetWarmstartSizeFilterSQP )
2177 {
2178  SCIPerrorMessage("method of filtersqp nonlinear solver is not implemented\n");
2179  SCIPABORT();
2180 
2181  return SCIP_OKAY; /*lint !e527*/
2182 } /*lint !e715*/
2183 
2184 /** stores warmstart information in buffer
2185  *
2186  * required size of buffer should have been obtained by SCIPnlpiGetWarmstartSize before
2187  *
2188  * input:
2189  * - nlpi datastructure for solver interface
2190  * - problem datastructure for problem instance
2191  * - buffer memory to store warmstart information
2192  *
2193  * output:
2194  * - buffer warmstart information in solver specific data structure
2195  */
2196 static
2197 SCIP_DECL_NLPIGETWARMSTARTMEMO( nlpiGetWarmstartMemoFilterSQP )
2198 {
2199  SCIPerrorMessage("method of filtersqp nonlinear solver is not implemented\n");
2200  SCIPABORT();
2201 
2202  return SCIP_OKAY; /*lint !e527*/
2203 } /*lint !e715*/
2204 
2205 /** sets warmstart information in solver
2206  *
2207  * write warmstart to buffer
2208  *
2209  * input:
2210  * - nlpi datastructure for solver interface
2211  * - problem datastructure for problem instance
2212  * - buffer warmstart information
2213  */
2214 static
2215 SCIP_DECL_NLPISETWARMSTARTMEMO( nlpiSetWarmstartMemoFilterSQP )
2216 {
2217  SCIPerrorMessage("method of filtersqp nonlinear solver is not implemented\n");
2218  SCIPABORT();
2219 
2220  return SCIP_OKAY; /*lint !e527*/
2221 } /*lint !e715*/
2222 
2223 /** gets integer parameter of NLP
2224  *
2225  * input:
2226  * - nlpi NLP interface structure
2227  * - problem datastructure for problem instance
2228  * - type parameter number
2229  * - ival pointer to store the parameter value
2230  *
2231  * output:
2232  * - ival parameter value
2233  */
2234 static
2235 SCIP_DECL_NLPIGETINTPAR( nlpiGetIntParFilterSQP )
2236 {
2237  assert(nlpi != NULL);
2238  assert(ival != NULL);
2239  assert(problem != NULL);
2240 
2241  switch( type )
2242  {
2244  {
2245  *ival = problem->fromscratch ? 1 : 0;
2246  break;
2247  }
2248 
2249  case SCIP_NLPPAR_VERBLEVEL:
2250  {
2251  *ival = problem->iprint;
2252  break;
2253  }
2254 
2255  case SCIP_NLPPAR_FEASTOL:
2256  {
2257  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
2258  return SCIP_PARAMETERWRONGTYPE;
2259  }
2260 
2261  case SCIP_NLPPAR_RELOBJTOL:
2262  {
2263  SCIPerrorMessage("relative objective tolerance parameter is of type real.\n");
2264  return SCIP_PARAMETERWRONGTYPE;
2265  }
2266 
2267  case SCIP_NLPPAR_LOBJLIM:
2268  {
2269  SCIPerrorMessage("objective limit parameter is of type real.\n");
2270  return SCIP_PARAMETERWRONGTYPE;
2271  }
2272 
2273  case SCIP_NLPPAR_INFINITY:
2274  {
2275  SCIPerrorMessage("infinity parameter is of type real.\n");
2276  return SCIP_PARAMETERWRONGTYPE;
2277  }
2278 
2279  case SCIP_NLPPAR_ITLIM:
2280  {
2281  *ival = problem->maxiter;
2282  break;
2283  }
2284 
2285  case SCIP_NLPPAR_TILIM:
2286  {
2287  SCIPerrorMessage("time limit parameter is of type real.\n");
2288  return SCIP_PARAMETERWRONGTYPE;
2289  }
2290 
2291  case SCIP_NLPPAR_OPTFILE:
2292  {
2293  SCIPerrorMessage("optfile parameter is of type string.\n");
2294  return SCIP_PARAMETERWRONGTYPE;
2295  }
2296 
2297  case SCIP_NLPPAR_FASTFAIL:
2298  {
2299  *ival = 0;
2300  break;
2301  }
2302 
2303  default:
2304  {
2305  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
2306  return SCIP_PARAMETERUNKNOWN;
2307  }
2308  }
2309 
2310  return SCIP_OKAY;
2311 } /*lint !e715*/
2312 
2313 /** sets integer parameter of NLP
2314  *
2315  * input:
2316  * - nlpi NLP interface structure
2317  * - problem datastructure for problem instance
2318  * - type parameter number
2319  * - ival parameter value
2320  */
2321 static
2322 SCIP_DECL_NLPISETINTPAR( nlpiSetIntParFilterSQP )
2323 {
2324  assert(nlpi != NULL);
2325  assert(problem != NULL);
2326 
2327  switch( type )
2328  {
2330  {
2331  if( ival == 0 || ival == 1 )
2332  {
2333  problem->fromscratch = (SCIP_Bool)ival;
2334  }
2335  else
2336  {
2337  SCIPerrorMessage("Value %d for parameter from scratch out of range {0, 1}\n", ival);
2338  return SCIP_PARAMETERWRONGVAL;
2339  }
2340  break;
2341  }
2342 
2343  case SCIP_NLPPAR_VERBLEVEL:
2344  {
2345  if( ival >= 0 )
2346  {
2347  problem->iprint = ival;
2348  }
2349  else
2350  {
2351  SCIPerrorMessage("Value %d for parameter from verbosity level out of range\n", ival);
2352  return SCIP_PARAMETERWRONGVAL;
2353  }
2354  break;
2355  }
2356 
2357  case SCIP_NLPPAR_FEASTOL:
2358  {
2359  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
2360  return SCIP_PARAMETERWRONGTYPE;
2361  }
2362 
2363  case SCIP_NLPPAR_RELOBJTOL:
2364  {
2365  SCIPerrorMessage("relative objective tolerance parameter is of type real.\n");
2366  return SCIP_PARAMETERWRONGTYPE;
2367  }
2368 
2369  case SCIP_NLPPAR_LOBJLIM:
2370  {
2371  SCIPerrorMessage("objective limit parameter is of type real.\n");
2372  return SCIP_PARAMETERWRONGTYPE;
2373  }
2374 
2375  case SCIP_NLPPAR_INFINITY:
2376  {
2377  SCIPerrorMessage("infinity parameter is of type real.\n");
2378  return SCIP_PARAMETERWRONGTYPE;
2379  }
2380 
2381  case SCIP_NLPPAR_ITLIM:
2382  {
2383  if( ival >= 0 )
2384  {
2385  problem->maxiter = ival;
2386  }
2387  else
2388  {
2389  SCIPerrorMessage("Value %d for parameter iteration limit is negative\n", ival);
2390  return SCIP_PARAMETERWRONGVAL;
2391  }
2392  break;
2393  }
2394 
2395  case SCIP_NLPPAR_TILIM:
2396  {
2397  SCIPerrorMessage("time limit parameter is of type real.\n");
2398  return SCIP_PARAMETERWRONGTYPE;
2399  }
2400 
2401  case SCIP_NLPPAR_OPTFILE:
2402  {
2403  SCIPerrorMessage("optfile parameter is of type string.\n");
2404  return SCIP_PARAMETERWRONGTYPE;
2405  }
2406 
2407  case SCIP_NLPPAR_FASTFAIL:
2408  {
2409  if( ival == 0 || ival == 1 )
2410  {
2411  SCIPdebugMessage("fast fail parameter not supported by FilterSQP interface yet. Ignored.\n");
2412  }
2413  else
2414  {
2415  SCIPerrorMessage("Value %d for parameter fastfail out of range {0, 1}\n", ival);
2416  return SCIP_PARAMETERWRONGVAL;
2417  }
2418  break;
2419  }
2420 
2421  default:
2422  {
2423  SCIPerrorMessage("Parameter %d not known to FilterSQP interface.\n", type);
2424  return SCIP_PARAMETERUNKNOWN;
2425  }
2426  }
2427 
2428  return SCIP_OKAY;
2429 } /*lint !e715*/
2430 
2431 /** gets floating point parameter of NLP
2432  *
2433  * input:
2434  * - nlpi NLP interface structure
2435  * - problem datastructure for problem instance, can be NULL only if type == SCIP_NLPPAR_INFINITY
2436  * - type parameter number
2437  * - dval pointer to store the parameter value
2438  *
2439  * output:
2440  * - dval parameter value
2441  */
2442 static
2443 SCIP_DECL_NLPIGETREALPAR( nlpiGetRealParFilterSQP )
2444 {
2445  assert(nlpi != NULL);
2446  assert(dval != NULL);
2447  assert(type == SCIP_NLPPAR_INFINITY || problem != NULL);
2448 
2449  switch( type )
2450  {
2452  {
2453  SCIPerrorMessage("from scratch parameter is of type int.\n");
2454  return SCIP_PARAMETERWRONGTYPE;
2455  }
2456 
2457  case SCIP_NLPPAR_VERBLEVEL:
2458  {
2459  SCIPerrorMessage("verbosity level parameter is of type int.\n");
2460  return SCIP_PARAMETERWRONGTYPE;
2461  }
2462 
2463  case SCIP_NLPPAR_FEASTOL:
2464  {
2465  *dval = problem->feastol;
2466  break;
2467  }
2468 
2469  case SCIP_NLPPAR_RELOBJTOL:
2470  {
2471  *dval = problem->opttol;
2472  break;
2473  }
2474 
2475  case SCIP_NLPPAR_LOBJLIM:
2476  {
2477  *dval = problem->fmin;
2478  break;
2479  }
2480 
2481  case SCIP_NLPPAR_INFINITY:
2482  {
2483  if( problem )
2484  {
2485  *dval = SCIPnlpiOracleGetInfinity(problem->oracle);
2486  }
2487  else
2488  {
2489  SCIP_NLPIDATA* data = SCIPnlpiGetData(nlpi);
2490  assert(data != NULL);
2491  *dval = data->infinity;
2492  }
2493  break;
2494  }
2495 
2496  case SCIP_NLPPAR_ITLIM:
2497  {
2498  SCIPerrorMessage("iteration limit parameter is of type int.\n");
2499  return SCIP_PARAMETERWRONGTYPE;
2500  }
2501 
2502  case SCIP_NLPPAR_TILIM:
2503  {
2504  *dval = problem->maxtime;
2505  break;
2506  }
2507 
2508  case SCIP_NLPPAR_OPTFILE:
2509  {
2510  SCIPerrorMessage("option file parameter is of type string.\n");
2511  return SCIP_PARAMETERWRONGTYPE;
2512  }
2513 
2514  case SCIP_NLPPAR_FASTFAIL:
2515  {
2516  SCIPerrorMessage("fastfail parameter is of type int.\n");
2517  return SCIP_PARAMETERWRONGTYPE;
2518  }
2519 
2520  default:
2521  {
2522  SCIPerrorMessage("Parameter %d not known to FilterSQP interface.\n", type);
2523  return SCIP_PARAMETERUNKNOWN;
2524  }
2525  }
2526 
2527  return SCIP_OKAY;
2528 } /*lint !e715*/
2529 
2530 /** sets floating point parameter of NLP
2531  *
2532  * input:
2533  * - nlpi NLP interface structure
2534  * - problem datastructure for problem instance, can be NULL only if type == SCIP_NLPPAR_INFINITY
2535  * - type parameter number
2536  * - dval parameter value
2537  */
2538 static
2539 SCIP_DECL_NLPISETREALPAR( nlpiSetRealParFilterSQP )
2540 {
2541  assert(nlpi != NULL);
2542  assert(type == SCIP_NLPPAR_INFINITY || problem != NULL);
2543 
2544  switch( type )
2545  {
2547  {
2548  SCIPerrorMessage("from scratch parameter is of type int.\n");
2549  return SCIP_PARAMETERWRONGTYPE;
2550  }
2551 
2552  case SCIP_NLPPAR_VERBLEVEL:
2553  {
2554  SCIPerrorMessage("verbosity level parameter is of type int.\n");
2555  return SCIP_PARAMETERWRONGTYPE;
2556  }
2557 
2558  case SCIP_NLPPAR_FEASTOL:
2559  {
2560  if( dval >= 0 )
2561  {
2562  problem->feastol = dval;
2563  }
2564  else
2565  {
2566  SCIPerrorMessage("Value %g for parameter feasibility tolerance is negative\n", dval);
2567  return SCIP_PARAMETERWRONGVAL;
2568  }
2569  break;
2570  }
2571 
2572  case SCIP_NLPPAR_RELOBJTOL:
2573  {
2574  if( dval >= 0 )
2575  {
2576  problem->opttol = dval;
2577  }
2578  else
2579  {
2580  SCIPerrorMessage("Value %g for parameter relative objective tolerance is negative\n", dval);
2581  return SCIP_PARAMETERWRONGVAL;
2582  }
2583  break;
2584  }
2585 
2586  case SCIP_NLPPAR_LOBJLIM:
2587  {
2588  problem->fmin = dval;
2589  break;
2590  }
2591 
2592  case SCIP_NLPPAR_INFINITY:
2593  {
2594  if( dval < 0.0 )
2595  return SCIP_PARAMETERWRONGVAL;
2596  if( problem )
2597  {
2598  SCIP_CALL( SCIPnlpiOracleSetInfinity(problem->oracle, dval) );
2599  }
2600  else
2601  {
2602  SCIP_NLPIDATA* data = SCIPnlpiGetData(nlpi);
2603  assert(data != NULL);
2604  data->infinity = dval;
2605  }
2606  break;
2607  }
2608 
2609  case SCIP_NLPPAR_ITLIM:
2610  {
2611  SCIPerrorMessage("iteration limit parameter is of type int.\n");
2612  return SCIP_PARAMETERWRONGTYPE;
2613  }
2614 
2615  case SCIP_NLPPAR_TILIM:
2616  {
2617  if( dval >= 0 )
2618  {
2619  problem->maxtime = dval;
2620  }
2621  else
2622  {
2623  SCIPerrorMessage("Value %g for parameter time limit is negative\n", dval);
2624  return SCIP_PARAMETERWRONGVAL;
2625  }
2626  break;
2627  }
2628 
2629  case SCIP_NLPPAR_OPTFILE:
2630  {
2631  SCIPerrorMessage("option file parameter is of type string.\n");
2632  return SCIP_PARAMETERWRONGTYPE;
2633  }
2634 
2635  case SCIP_NLPPAR_FASTFAIL:
2636  {
2637  SCIPerrorMessage("fastfail parameter is of type int.\n");
2638  return SCIP_PARAMETERWRONGTYPE;
2639  }
2640 
2641  default:
2642  {
2643  SCIPerrorMessage("Parameter %d not known to FilterSQP interface.\n", type);
2644  return SCIP_PARAMETERUNKNOWN;
2645  }
2646  }
2647 
2648  return SCIP_OKAY;
2649 } /*lint !e715*/
2650 
2651 /** gets string parameter of NLP
2652  *
2653  * input:
2654  * - nlpi NLP interface structure
2655  * - problem datastructure for problem instance
2656  * - type parameter number
2657  * - sval pointer to store the string value, the user must not modify the string
2658  *
2659  * output:
2660  * - sval parameter value
2661  */
2662 static
2663 SCIP_DECL_NLPIGETSTRINGPAR( nlpiGetStringParFilterSQP )
2664 {
2665  assert(nlpi != NULL);
2666  assert(problem != NULL);
2667 
2668  switch( type )
2669  {
2671  {
2672  SCIPerrorMessage("from scratch parameter is of type int.\n");
2673  return SCIP_PARAMETERWRONGTYPE;
2674  }
2675 
2676  case SCIP_NLPPAR_VERBLEVEL:
2677  {
2678  SCIPerrorMessage("verbosity level parameter is of type int.\n");
2679  return SCIP_PARAMETERWRONGTYPE;
2680  }
2681 
2682  case SCIP_NLPPAR_FEASTOL:
2683  {
2684  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
2685  return SCIP_PARAMETERWRONGTYPE;
2686  }
2687 
2688  case SCIP_NLPPAR_RELOBJTOL:
2689  {
2690  SCIPerrorMessage("objective tolerance parameter is of type real.\n");
2691  return SCIP_PARAMETERWRONGTYPE;
2692  }
2693 
2694  case SCIP_NLPPAR_LOBJLIM:
2695  {
2696  SCIPerrorMessage("objective limit parameter is of type real.\n");
2697  return SCIP_PARAMETERWRONGTYPE;
2698  }
2699 
2700  case SCIP_NLPPAR_INFINITY:
2701  {
2702  SCIPerrorMessage("infinity parameter is of type real.\n");
2703  return SCIP_PARAMETERWRONGTYPE;
2704  }
2705 
2706  case SCIP_NLPPAR_ITLIM:
2707  {
2708  SCIPerrorMessage("iteration limit parameter is of type int.\n");
2709  return SCIP_PARAMETERWRONGTYPE;
2710  }
2711 
2712  case SCIP_NLPPAR_TILIM:
2713  {
2714  SCIPerrorMessage("time limit parameter is of type real.\n");
2715  return SCIP_PARAMETERWRONGTYPE;
2716  }
2717 
2718  case SCIP_NLPPAR_OPTFILE:
2719  {
2720  *sval = NULL;
2721  return SCIP_OKAY;
2722  }
2723 
2724  case SCIP_NLPPAR_FASTFAIL:
2725  {
2726  SCIPerrorMessage("fastfail parameter is of type int.\n");
2727  return SCIP_PARAMETERWRONGTYPE;
2728  }
2729 
2730  default:
2731  {
2732  SCIPerrorMessage("Parameter %d not known to FilterSQP interface.\n", type);
2733  return SCIP_PARAMETERUNKNOWN;
2734  }
2735  }
2736 } /*lint !e715*/
2737 
2738 /** sets string parameter of NLP
2739  *
2740  * input:
2741  * - nlpi NLP interface structure
2742  * - problem datastructure for problem instance
2743  * - type parameter number
2744  * - sval parameter value
2745  */
2746 static
2747 SCIP_DECL_NLPISETSTRINGPAR( nlpiSetStringParFilterSQP )
2748 {
2749  assert(nlpi != NULL);
2750  assert(problem != NULL);
2751 
2752  switch( type )
2753  {
2755  {
2756  SCIPerrorMessage("from scratch parameter is of type int.\n");
2757  return SCIP_PARAMETERWRONGTYPE;
2758  }
2759 
2760  case SCIP_NLPPAR_VERBLEVEL:
2761  {
2762  SCIPerrorMessage("verbosity level parameter is of type int.\n");
2763  return SCIP_PARAMETERWRONGTYPE;
2764  }
2765 
2766  case SCIP_NLPPAR_FEASTOL:
2767  {
2768  SCIPerrorMessage("feasibility tolerance parameter is of type real.\n");
2769  return SCIP_PARAMETERWRONGTYPE;
2770  }
2771 
2772  case SCIP_NLPPAR_RELOBJTOL:
2773  {
2774  SCIPerrorMessage("objective tolerance parameter is of type real.\n");
2775  return SCIP_PARAMETERWRONGTYPE;
2776  }
2777 
2778  case SCIP_NLPPAR_LOBJLIM:
2779  {
2780  SCIPerrorMessage("objective limit parameter is of type real.\n");
2781  return SCIP_PARAMETERWRONGTYPE;
2782  }
2783 
2784  case SCIP_NLPPAR_INFINITY:
2785  {
2786  SCIPerrorMessage("infinity parameter is of type real.\n");
2787  return SCIP_PARAMETERWRONGTYPE;
2788  }
2789 
2790  case SCIP_NLPPAR_ITLIM:
2791  {
2792  SCIPerrorMessage("iteration limit parameter is of type int.\n");
2793  return SCIP_PARAMETERWRONGTYPE;
2794  }
2795 
2796  case SCIP_NLPPAR_TILIM:
2797  {
2798  SCIPerrorMessage("time limit parameter is of type real.\n");
2799  return SCIP_PARAMETERWRONGTYPE;
2800  }
2801 
2802  case SCIP_NLPPAR_OPTFILE:
2803  {
2804  SCIPmessagePrintWarning(SCIPnlpiGetData(nlpi)->messagehdlr, "Parameter optfile not supported by FilterSQP interface. Ignored.\n");
2805  return SCIP_OKAY;
2806  }
2807 
2808  case SCIP_NLPPAR_FASTFAIL:
2809  {
2810  SCIPerrorMessage("fastfail parameter is of type int.\n");
2811  return SCIP_PARAMETERWRONGTYPE;
2812  }
2813 
2814  default:
2815  {
2816  SCIPerrorMessage("Parameter %d not known to Ipopt interface.\n", type);
2817  return SCIP_PARAMETERUNKNOWN;
2818  }
2819  }
2820 } /*lint !e715*/
2821 
2822 /** sets message handler for message output
2823  *
2824  * input:
2825  * - nlpi NLP interface structure
2826  * - messagehdlr SCIP message handler, or NULL to suppress all output
2827  */
2828 static
2829 SCIP_DECL_NLPISETMESSAGEHDLR( nlpiSetMessageHdlrFilterSQP )
2830 {
2831  SCIP_NLPIDATA* nlpidata;
2832 
2833  assert(nlpi != NULL);
2834 
2835  nlpidata = SCIPnlpiGetData(nlpi);
2836  assert(nlpidata != NULL);
2837 
2838  nlpidata->messagehdlr = messagehdlr;
2839 
2840  return SCIP_OKAY; /*lint !e527*/
2841 } /*lint !e715*/
2842 
2843 /*
2844  * NLP solver interface specific interface methods
2845  */
2846 
2847 /** create solver interface for FilterSQP solver */
2849  BMS_BLKMEM* blkmem, /**< block memory data structure */
2850  SCIP_NLPI** nlpi /**< pointer to buffer for nlpi address */
2851  )
2852 {
2853  SCIP_NLPIDATA* nlpidata;
2854 
2855  assert(blkmem != NULL);
2856  assert(nlpi != NULL);
2857 
2858  /* create filterSQP solver interface data */
2859  SCIP_ALLOC( BMSallocBlockMemory(blkmem, &nlpidata) );
2860 
2861  nlpidata->blkmem = blkmem;
2862  nlpidata->messagehdlr = NULL;
2863  nlpidata->infinity = SCIP_DEFAULT_INFINITY;
2864  nlpidata->randnumgen = NULL;
2865 
2866  /* create solver interface */
2867  SCIP_CALL( SCIPnlpiCreate(nlpi,
2869  nlpiCopyFilterSQP, nlpiFreeFilterSQP, nlpiGetSolverPointerFilterSQP,
2870  nlpiCreateProblemFilterSQP, nlpiFreeProblemFilterSQP, nlpiGetProblemPointerFilterSQP,
2871  nlpiAddVarsFilterSQP, nlpiAddConstraintsFilterSQP, nlpiSetObjectiveFilterSQP,
2872  nlpiChgVarBoundsFilterSQP, nlpiChgConsSidesFilterSQP, nlpiDelVarSetFilterSQP, nlpiDelConstraintSetFilterSQP,
2873  nlpiChgLinearCoefsFilterSQP, nlpiChgQuadraticCoefsFilterSQP, nlpiChgExprtreeFilterSQP, nlpiChgNonlinCoefFilterSQP,
2874  nlpiChgObjConstantFilterSQP, nlpiSetInitialGuessFilterSQP, nlpiSolveFilterSQP, nlpiGetSolstatFilterSQP, nlpiGetTermstatFilterSQP,
2875  nlpiGetSolutionFilterSQP, nlpiGetStatisticsFilterSQP,
2876  nlpiGetWarmstartSizeFilterSQP, nlpiGetWarmstartMemoFilterSQP, nlpiSetWarmstartMemoFilterSQP,
2877  nlpiGetIntParFilterSQP, nlpiSetIntParFilterSQP, nlpiGetRealParFilterSQP, nlpiSetRealParFilterSQP, nlpiGetStringParFilterSQP, nlpiSetStringParFilterSQP,
2878  nlpiSetMessageHdlrFilterSQP,
2879  nlpidata) );
2880 
2881  return SCIP_OKAY;
2882 }
2883 
2884 /** gets string that identifies filterSQP (version number) */
2886  void
2887  )
2888 {
2889  return "filterSQP"; /* TODO version number? */
2890 }
2891 
2892 /** gets string that describes filterSQP */
2894  void
2895  )
2896 {
2897  return NLPI_DESC;
2898 }
2899 
2900 /** returns whether filterSQP is available, i.e., whether it has been linked in */
2902  void
2903  )
2904 {
2905  return TRUE;
2906 }
static SCIP_DECL_NLPIFREEPROBLEM(nlpiFreeProblemFilterSQP)
static SCIP_DECL_NLPIGETWARMSTARTSIZE(nlpiGetWarmstartSizeFilterSQP)
#define BMSfreeBlockMemoryArrayNull(mem, ptr, num)
Definition: memory.h:449
enum SCIP_NlpTermStat SCIP_NLPTERMSTAT
Definition: type_nlpi.h:84
static SCIP_DECL_NLPIGETSOLSTAT(nlpiGetSolstatFilterSQP)
const SCIP_Real * SCIPnlpiOracleGetVarUbs(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2225
const char * SCIPgetSolverNameFilterSQP(void)
SCIP_RETCODE SCIPnlpiOracleGetJacobianSparsity(SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
Definition: nlpioracle.c:2510
static SCIP_RETCODE setupGradients(BMS_BLKMEM *blkmem, SCIP_NLPIORACLE *oracle, fint **la, int *lasize, real **a)
SCIP_RETCODE SCIPcreateNlpSolverFilterSQP(BMS_BLKMEM *blkmem, SCIP_NLPI **nlpi)
SCIP_RETCODE SCIPnlpiOracleGetHessianLagSparsity(SCIP_NLPIORACLE *oracle, const int **offset, const int **col)
Definition: nlpioracle.c:2790
int SCIPnlpiOracleGetNVars(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2195
static SCIP_DECL_NLPIGETWARMSTARTMEMO(nlpiGetWarmstartMemoFilterSQP)
#define infinity
Definition: gastrans.c:71
static SCIP_DECL_NLPIADDCONSTRAINTS(nlpiAddConstraintsFilterSQP)
static SCIP_DECL_NLPICHGVARBOUNDS(nlpiChgVarBoundsFilterSQP)
static SCIP_DECL_NLPISETWARMSTARTMEMO(nlpiSetWarmstartMemoFilterSQP)
#define NLPI_PRIORITY
fint n_bqpd_calls
internal methods for NLPI solver interfaces
SCIP_Bool SCIPisFilterSQPAvailableFilterSQP(void)
SCIP_RETCODE SCIPnlpiOracleSetObjective(SCIP_NLPIORACLE *oracle, const SCIP_Real constant, int nlin, const int *lininds, const SCIP_Real *linvals, int nquadelems, const SCIP_QUADELEM *quadelems, const int *exprvaridxs, const SCIP_EXPRTREE *exprtree)
Definition: nlpioracle.c:1607
SCIP_Real * varubdualvalues
static SCIP_DECL_NLPISOLVE(nlpiSolveFilterSQP)
methods to store an NLP and request function, gradient, and hessian values
SCIP_NLPIORACLE * oracle
#define FALSE
Definition: def.h:64
SCIP_Real SCIPnlpiOracleGetConstraintLhs(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2276
static SCIP_DECL_NLPICHGCONSSIDES(nlpiChgConsSidesFilterSQP)
#define TRUE
Definition: def.h:63
#define NLPI_NAME
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:53
void F77_FUNC(filtersqp, FILTERSQP)
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveValue(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Real *objval)
Definition: nlpioracle.c:2397
static SCIP_DECL_NLPIADDVARS(nlpiAddVarsFilterSQP)
SCIP_Real * primalvalues
#define WORKSPACEGROWTHFACTOR
static SCIP_TIME gettime(void)
#define BMSallocMemoryArray(ptr, num)
Definition: memory.h:105
static SCIP_DECL_NLPIGETSOLUTION(nlpiGetSolutionFilterSQP)
#define SCIPdebugMessage
Definition: pub_message.h:77
SCIP_RETCODE SCIPnlpiSetMessageHdlr(SCIP_NLPI *nlpi, SCIP_MESSAGEHDLR *messagehdlr)
Definition: nlpi.c:720
#define OPTTOLFACTOR
static SCIP_DECL_NLPIGETPROBLEMPOINTER(nlpiGetProblemPointerFilterSQP)
static pthread_mutex_t filtersqpmutex
SCIP_RETCODE SCIPnlpiOracleEvalObjectiveGradient(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *objval, SCIP_Real *objgrad)
Definition: nlpioracle.c:2461
SCIP_RETCODE SCIPnlpiOracleDelConsSet(SCIP_NLPIORACLE *oracle, int *delstats)
Definition: nlpioracle.c:1824
SCIP_RETCODE SCIPnlpiOracleChgObjConstant(SCIP_NLPIORACLE *oracle, SCIP_Real objconstant)
Definition: nlpioracle.c:2179
filterSQP NLP interface
int SCIPnlpiOracleGetNConstraints(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2205
real eps
SCIP_Real * varlbdualvalues
static SCIP_DECL_NLPIDELVARSET(nlpiDelVarSetFilterSQP)
static SCIP_DECL_NLPIGETTERMSTAT(nlpiGetTermstatFilterSQP)
long ftnlen
static SCIP_DECL_NLPISETINITIALGUESS(nlpiSetInitialGuessFilterSQP)
int SCIPnlpiOracleGetConstraintDegree(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2317
#define BMSfreeMemoryArray(ptr)
Definition: memory.h:129
SCIP_RETCODE SCIPnlpiOracleChgExprtree(SCIP_NLPIORACLE *oracle, int considx, const int *exprvaridxs, const SCIP_EXPRTREE *exprtree)
Definition: nlpioracle.c:2097
static SCIP_DECL_NLPICHGQUADCOEFS(nlpiChgQuadraticCoefsFilterSQP)
SCIP_RETCODE SCIPnlpiOracleCreate(BMS_BLKMEM *blkmem, SCIP_NLPIORACLE **oracle)
Definition: nlpioracle.c:1325
#define DEFAULT_MAXITER
#define SCIPerrorMessage
Definition: pub_message.h:45
static SCIP_DECL_NLPIGETSOLVERPOINTER(nlpiGetSolverPointerFilterSQP)
const SCIP_Real * SCIPnlpiOracleGetVarLbs(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:2215
SCIP_RETCODE SCIPnlpiOracleChgVarBounds(SCIP_NLPIORACLE *oracle, int nvars, const int *indices, const SCIP_Real *lbs, const SCIP_Real *ubs)
Definition: nlpioracle.c:1644
struct SCIP_NlpiData SCIP_NLPIDATA
Definition: type_nlpi.h:38
enum SCIP_NlpSolStat SCIP_NLPSOLSTAT
Definition: type_nlpi.h:69
SCIP_RETCODE SCIPnlpiOracleAddVars(SCIP_NLPIORACLE *oracle, int nvars, const SCIP_Real *lbs, const SCIP_Real *ubs, const char **varnames)
Definition: nlpioracle.c:1444
SCIP_Bool warmstart
fint scale_mode
void SCIPmessagePrintWarning(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr,...)
Definition: message.c:417
#define BMSallocClearBlockMemoryArray(mem, ptr, num)
Definition: memory.h:436
internal miscellaneous methods
int fint
void SCIPrandomFree(SCIP_RANDNUMGEN **randnumgen, BMS_BLKMEM *blkmem)
Definition: misc.c:9350
fint n_bqpd_prfint
static SCIP_DECL_NLPIGETINTPAR(nlpiGetIntParFilterSQP)
#define SCIP_CALL(x)
Definition: def.h:350
static SCIP_DECL_NLPISETINTPAR(nlpiSetIntParFilterSQP)
static SCIP_DECL_NLPICOPY(nlpiCopyFilterSQP)
const char * SCIPgetSolverDescFilterSQP(void)
void SCIPmessagePrintInfo(SCIP_MESSAGEHDLR *messagehdlr, const char *formatstr,...)
Definition: message.c:584
#define MAXNRUNS
static SCIP_DECL_NLPIFREE(nlpiFreeFilterSQP)
static SCIP_DECL_NLPIGETREALPAR(nlpiGetRealParFilterSQP)
SCIP_Real SCIPnlpiOracleGetConstraintRhs(SCIP_NLPIORACLE *oracle, int considx)
Definition: nlpioracle.c:2289
SCIP_Real SCIPnlpiOracleGetInfinity(SCIP_NLPIORACLE *oracle)
Definition: nlpioracle.c:1397
SCIP_Real solvetime
SCIP_RETCODE SCIPnlpiOracleFree(SCIP_NLPIORACLE **oracle)
Definition: nlpioracle.c:1352
static SCIP_Real timeelapsed(SCIP_NLPIDATA *nlpidata)
real infty
#define BMSfreeBlockMemory(mem, ptr)
Definition: memory.h:446
fint phr
SCIP_RETCODE SCIPnlpiCreate(SCIP_NLPI **nlpi, const char *name, const char *description, int priority, SCIP_DECL_NLPICOPY((*nlpicopy)), SCIP_DECL_NLPIFREE((*nlpifree)), SCIP_DECL_NLPIGETSOLVERPOINTER((*nlpigetsolverpointer)), SCIP_DECL_NLPICREATEPROBLEM((*nlpicreateproblem)), SCIP_DECL_NLPIFREEPROBLEM((*nlpifreeproblem)), SCIP_DECL_NLPIGETPROBLEMPOINTER((*nlpigetproblempointer)), SCIP_DECL_NLPIADDVARS((*nlpiaddvars)), SCIP_DECL_NLPIADDCONSTRAINTS((*nlpiaddconstraints)), SCIP_DECL_NLPISETOBJECTIVE((*nlpisetobjective)), SCIP_DECL_NLPICHGVARBOUNDS((*nlpichgvarbounds)), SCIP_DECL_NLPICHGCONSSIDES((*nlpichgconssides)), SCIP_DECL_NLPIDELVARSET((*nlpidelvarset)), SCIP_DECL_NLPIDELCONSSET((*nlpidelconsset)), SCIP_DECL_NLPICHGLINEARCOEFS((*nlpichglinearcoefs)), SCIP_DECL_NLPICHGQUADCOEFS((*nlpichgquadcoefs)), SCIP_DECL_NLPICHGEXPRTREE((*nlpichgexprtree)), SCIP_DECL_NLPICHGNONLINCOEF((*nlpichgnonlincoef)), SCIP_DECL_NLPICHGOBJCONSTANT((*nlpichgobjconstant)), SCIP_DECL_NLPISETINITIALGUESS((*nlpisetinitialguess)), SCIP_DECL_NLPISOLVE((*nlpisolve)), SCIP_DECL_NLPIGETSOLSTAT((*nlpigetsolstat)), SCIP_DECL_NLPIGETTERMSTAT((*nlpigettermstat)), SCIP_DECL_NLPIGETSOLUTION((*nlpigetsolution)), SCIP_DECL_NLPIGETSTATISTICS((*nlpigetstatistics)), SCIP_DECL_NLPIGETWARMSTARTSIZE((*nlpigetwarmstartsize)), SCIP_DECL_NLPIGETWARMSTARTMEMO((*nlpigetwarmstartmemo)), SCIP_DECL_NLPISETWARMSTARTMEMO((*nlpisetwarmstartmemo)), SCIP_DECL_NLPIGETINTPAR((*nlpigetintpar)), SCIP_DECL_NLPISETINTPAR((*nlpisetintpar)), SCIP_DECL_NLPIGETREALPAR((*nlpigetrealpar)), SCIP_DECL_NLPISETREALPAR((*nlpisetrealpar)), SCIP_DECL_NLPIGETSTRINGPAR((*nlpigetstringpar)), SCIP_DECL_NLPISETSTRINGPAR((*nlpisetstringpar)), SCIP_DECL_NLPISETMESSAGEHDLR((*nlpisetmessagehdlr)), SCIP_NLPIDATA *nlpidata)
Definition: nlpi.c:40
#define MINEPS
static SCIP_Bool timelimitreached(SCIP_NLPIDATA *nlpidata, SCIP_NLPIPROBLEM *nlpiproblem)
#define SCIP_Bool
Definition: def.h:61
#define BMSallocBlockMemoryArray(mem, ptr, num)
Definition: memory.h:435
SCIP_RETCODE SCIPnlpiOracleSetProblemName(SCIP_NLPIORACLE *oracle, const char *name)
Definition: nlpioracle.c:1409
fint phl
#define RANDSEED
SCIP_RETCODE SCIPnlpiOracleAddConstraints(SCIP_NLPIORACLE *oracle, int nconss, const SCIP_Real *lhss, const SCIP_Real *rhss, const int *nlininds, int *const *lininds, SCIP_Real *const *linvals, const int *nquadelems, SCIP_QUADELEM *const *quadelems, int *const *exprvaridxs, SCIP_EXPRTREE *const *exprtrees, const char **consnames)
Definition: nlpioracle.c:1529
real ubd
#define BMSfreeBlockMemoryArray(mem, ptr, num)
Definition: memory.h:448
void SCIPnlpStatisticsSetNIterations(SCIP_NLPSTATISTICS *statistics, int niterations)
Definition: nlpi.c:836
#define MAX(x, y)
Definition: tclique_def.h:75
void SCIPnlpStatisticsSetTotalTime(SCIP_NLPSTATISTICS *statistics, SCIP_Real totaltime)
Definition: nlpi.c:846
static int calcGrowSize(int num)
SCIP_RETCODE SCIPnlpiOracleEvalConstraintValue(SCIP_NLPIORACLE *oracle, int considx, const SCIP_Real *x, SCIP_Real *conval)
Definition: nlpioracle.c:2416
static SCIP_DECL_NLPIGETSTRINGPAR(nlpiGetStringParFilterSQP)
SCIP_RETCODE SCIPnlpiOracleChgLinearCoefs(SCIP_NLPIORACLE *oracle, int considx, int nentries, const int *varidxs, const SCIP_Real *newcoefs)
Definition: nlpioracle.c:1902
#define BMScopyMemoryArray(ptr, source, num)
Definition: memory.h:116
SCIP_RETCODE SCIPnlpiOracleChgQuadCoefs(SCIP_NLPIORACLE *oracle, int considx, int nquadelems, const SCIP_QUADELEM *quadelems)
Definition: nlpioracle.c:1999
#define NLPI_DESC
static SCIP_DECL_NLPISETMESSAGEHDLR(nlpiSetMessageHdlrFilterSQP)
SCIP_Real * initguess
SCIP_RETCODE SCIPnlpiOracleChgConsSides(SCIP_NLPIORACLE *oracle, int nconss, const int *indices, const SCIP_Real *lhss, const SCIP_Real *rhss)
Definition: nlpioracle.c:1680
#define BMSclearMemory(ptr)
Definition: memory.h:111
#define DEFAULT_LOBJLIM
SCIP_Bool fromscratch
static SCIP_DECL_NLPIGETSTATISTICS(nlpiGetStatisticsFilterSQP)
static SCIP_DECL_NLPICHGNONLINCOEF(nlpiChgNonlinCoefFilterSQP)
#define DEFAULT_FEASOPTTOL
SCIP_Real SCIPrandomGetReal(SCIP_RANDNUMGEN *randnumgen, SCIP_Real minrandval, SCIP_Real maxrandval)
Definition: misc.c:9388
SCIP_RETCODE SCIPnlpiOracleDelVarSet(SCIP_NLPIORACLE *oracle, int *delstats)
Definition: nlpioracle.c:1714
#define SCIP_DEFAULT_INFINITY
Definition: def.h:154
static SCIP_RETCODE processSolveOutcome(SCIP_NLPIDATA *nlpidata, SCIP_NLPIPROBLEM *problem, fint ifail, real *x, real *lam)
static SCIP_DECL_NLPIDELCONSSET(nlpiDelConstraintSetFilterSQP)
time_t sec
SCIP_RETCODE SCIPnlpiOracleChgExprParam(SCIP_NLPIORACLE *oracle, int considx, int paramidx, SCIP_Real paramval)
Definition: nlpioracle.c:2154
SCIP_NLPIDATA * SCIPnlpiGetData(SCIP_NLPI *nlpi)
Definition: nlpi.c:733
static SCIP_DECL_NLPICHGOBJCONSTANT(nlpiChgObjConstantFilterSQP)
SCIP_RETCODE SCIPnlpiOracleEvalHessianLag(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real objfactor, const SCIP_Real *lambda, SCIP_Real *hessian)
Definition: nlpioracle.c:2889
#define MAXPERTURB
real tt
#define SCIP_Real
Definition: def.h:149
static SCIP_DECL_NLPICREATEPROBLEM(nlpiCreateProblemFilterSQP)
SCIP_RETCODE SCIPrandomCreate(SCIP_RANDNUMGEN **randnumgen, BMS_BLKMEM *blkmem, unsigned int initialseed)
Definition: misc.c:9334
static void invalidateSolution(SCIP_NLPIPROBLEM *problem)
static SCIP_RETCODE setupHessian(BMS_BLKMEM *blkmem, SCIP_NLPIORACLE *oracle, fint **la, int *lasize)
#define SCIP_INVALID
Definition: def.h:169
static SCIP_DECL_NLPICHGEXPRTREE(nlpiChgExprtreeFilterSQP)
static SCIP_RETCODE setupStart(SCIP_NLPIDATA *data, SCIP_NLPIPROBLEM *problem, real *x, SCIP_Bool *success)
SCIP_RETCODE SCIPnlpiOracleSetInfinity(SCIP_NLPIORACLE *oracle, SCIP_Real infinity)
Definition: nlpioracle.c:1381
#define SCIPisFinite(x)
Definition: pub_misc.h:1768
SCIP_Real * consdualvalues
SCIP_NLPTERMSTAT termstat
static SCIP_DECL_NLPISETREALPAR(nlpiSetRealParFilterSQP)
SCIP_RETCODE SCIPnlpiOracleEvalJacobian(SCIP_NLPIORACLE *oracle, const SCIP_Real *x, SCIP_Bool isnewx, SCIP_Real *convals, SCIP_Real *jacobi)
Definition: nlpioracle.c:2651
SCIP_NLPSOLSTAT solstat
#define BMSallocBlockMemory(mem, ptr)
Definition: memory.h:433
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:112
static SCIP_DECL_NLPISETSTRINGPAR(nlpiSetStringParFilterSQP)
struct BMS_BlkMem BMS_BLKMEM
Definition: memory.h:419
static SCIP_DECL_NLPISETOBJECTIVE(nlpiSetObjectiveFilterSQP)
#define SCIP_ALLOC(x)
Definition: def.h:361
#define SCIPABORT()
Definition: def.h:322
fint phc
fint phe
double real
#define BMSreallocBlockMemoryArray(mem, ptr, oldnum, newnum)
Definition: memory.h:439
static SCIP_DECL_NLPICHGLINEARCOEFS(nlpiChgLinearCoefsFilterSQP)
SCIP_RETCODE SCIPnlpiSetRealPar(SCIP_NLPI *nlpi, SCIP_NLPIPROBLEM *problem, SCIP_NLPPARAM type, SCIP_Real dval)
Definition: nlpi.c:671