1 // Copyright (C) 2004, 2010 International Business Machines and others.
2 // All Rights Reserved.
3 // This code is published under the Eclipse Public License.
4 //
5 // $Id$
6 //
7 // Authors:  Carl Laird, Andreas Waechter     IBM    2004-08-13
8 
9 #ifndef __IPORIGIPOPTNLP_HPP__
10 #define __IPORIGIPOPTNLP_HPP__
11 
12 #include "IpIpoptNLP.hpp"
13 #include "IpException.hpp"
14 #include "IpTimingStatistics.hpp"
15 
16 namespace Ipopt
17 {
18 
19   /** enumeration for the Hessian information type. */
20   enum HessianApproximationType {
21     EXACT=0,
22     LIMITED_MEMORY
23   };
24 
25   /** enumeration for the Hessian approximation space. */
26   enum HessianApproximationSpace {
27     NONLINEAR_VARS=0,
28     ALL_VARS
29   };
30 
31   /** This class maps the traditional NLP into
32    *  something that is more useful by Ipopt.
33    *  This class takes care of storing the
34    *  calculated model results, handles caching,
35    *  and (some day) takes care of addition of slacks.
36    */
37   class OrigIpoptNLP : public IpoptNLP
38   {
39   public:
40     /**@name Constructors/Destructors */
41     //@{
42     OrigIpoptNLP(const SmartPtr<const Journalist>& jnlst,
43                  const SmartPtr<NLP>& nlp,
44                  const SmartPtr<NLPScalingObject>& nlp_scaling);
45 
46     /** Default destructor */
47     virtual ~OrigIpoptNLP();
48     //@}
49 
50     /** Initialize - overloaded from IpoptNLP */
51     virtual bool Initialize(const Journalist& jnlst,
52                             const OptionsList& options,
53                             const std::string& prefix);
54 
55     /** Initialize (create) structures for
56      *  the iteration data */
57     virtual bool InitializeStructures(SmartPtr<Vector>& x,
58                                       bool init_x,
59                                       SmartPtr<Vector>& y_c,
60                                       bool init_y_c,
61                                       SmartPtr<Vector>& y_d,
62                                       bool init_y_d,
63                                       SmartPtr<Vector>& z_L,
64                                       bool init_z_L,
65                                       SmartPtr<Vector>& z_U,
66                                       bool init_z_U,
67                                       SmartPtr<Vector>& v_L,
68                                       SmartPtr<Vector>& v_U
69                                      );
70 
71     /** Method accessing the GetWarmStartIterate of the NLP */
GetWarmStartIterate(IteratesVector & warm_start_iterate)72     virtual bool GetWarmStartIterate(IteratesVector& warm_start_iterate)
73     {
74       return nlp_->GetWarmStartIterate(warm_start_iterate);
75     }
76     /** Accessor methods for model data */
77     //@{
78     /** Objective value */
79     virtual Number f(const Vector& x);
80 
81     /** Objective value (depending in mu) - incorrect version for
82      *  OrigIpoptNLP */
83     virtual Number f(const Vector& x, Number mu);
84 
85     /** Gradient of the objective */
86     virtual SmartPtr<const Vector> grad_f(const Vector& x);
87 
88     /** Gradient of the objective (depending in mu) - incorrect
89      *  version for OrigIpoptNLP */
90     virtual SmartPtr<const Vector> grad_f(const Vector& x, Number mu);
91 
92     /** Equality constraint residual */
93     virtual SmartPtr<const Vector> c(const Vector& x);
94 
95     /** Jacobian Matrix for equality constraints */
96     virtual SmartPtr<const Matrix> jac_c(const Vector& x);
97 
98     /** Inequality constraint residual (reformulated
99      *  as equalities with slacks */
100     virtual SmartPtr<const Vector> d(const Vector& x);
101 
102     /** Jacobian Matrix for inequality constraints*/
103     virtual SmartPtr<const Matrix> jac_d(const Vector& x);
104 
105     /** Hessian of the Lagrangian */
106     virtual SmartPtr<const SymMatrix> h(const Vector& x,
107                                         Number obj_factor,
108                                         const Vector& yc,
109                                         const Vector& yd
110                                        );
111 
112     /** Hessian of the Lagrangian (depending in mu) - incorrect
113      *  version for OrigIpoptNLP */
114     virtual SmartPtr<const SymMatrix> h(const Vector& x,
115                                         Number obj_factor,
116                                         const Vector& yc,
117                                         const Vector& yd,
118                                         Number mu);
119 
120     /** Provides a Hessian matrix from the correct matrix space with
121      *  uninitialized values.  This can be used in LeastSquareMults to
122      *  obtain a "zero Hessian". */
123     virtual SmartPtr<const SymMatrix> uninitialized_h();
124 
125     /** Lower bounds on x */
x_L() const126     virtual SmartPtr<const Vector> x_L() const
127     {
128       return x_L_;
129     }
130 
131     /** Permutation matrix (x_L_ -> x) */
Px_L() const132     virtual SmartPtr<const Matrix> Px_L() const
133     {
134       return Px_L_;
135     }
136 
137     /** Upper bounds on x */
x_U() const138     virtual SmartPtr<const Vector> x_U() const
139     {
140       return x_U_;
141     }
142 
143     /** Permutation matrix (x_U_ -> x */
Px_U() const144     virtual SmartPtr<const Matrix> Px_U() const
145     {
146       return Px_U_;
147     }
148 
149     /** Lower bounds on d */
d_L() const150     virtual SmartPtr<const Vector> d_L() const
151     {
152       return d_L_;
153     }
154 
155     /** Permutation matrix (d_L_ -> d) */
Pd_L() const156     virtual SmartPtr<const Matrix> Pd_L() const
157     {
158       return Pd_L_;
159     }
160 
161     /** Upper bounds on d */
d_U() const162     virtual SmartPtr<const Vector> d_U() const
163     {
164       return d_U_;
165     }
166 
167     /** Permutation matrix (d_U_ -> d */
Pd_U() const168     virtual SmartPtr<const Matrix> Pd_U() const
169     {
170       return Pd_U_;
171     }
172 
HessianMatrixSpace() const173     virtual SmartPtr<const SymMatrixSpace> HessianMatrixSpace() const
174     {
175       return h_space_;
176     }
177 
x_space() const178     virtual SmartPtr<const VectorSpace> x_space() const
179     {
180       return x_space_;
181     }
182     //@}
183 
184     /** Accessor method for vector/matrix spaces pointers */
185     virtual void GetSpaces(SmartPtr<const VectorSpace>& x_space,
186                            SmartPtr<const VectorSpace>& c_space,
187                            SmartPtr<const VectorSpace>& d_space,
188                            SmartPtr<const VectorSpace>& x_l_space,
189                            SmartPtr<const MatrixSpace>& px_l_space,
190                            SmartPtr<const VectorSpace>& x_u_space,
191                            SmartPtr<const MatrixSpace>& px_u_space,
192                            SmartPtr<const VectorSpace>& d_l_space,
193                            SmartPtr<const MatrixSpace>& pd_l_space,
194                            SmartPtr<const VectorSpace>& d_u_space,
195                            SmartPtr<const MatrixSpace>& pd_u_space,
196                            SmartPtr<const MatrixSpace>& Jac_c_space,
197                            SmartPtr<const MatrixSpace>& Jac_d_space,
198                            SmartPtr<const SymMatrixSpace>& Hess_lagrangian_space);
199 
200     /** Method for adapting the variable bounds.  This is called if
201      *  slacks are becoming too small */
202     virtual void AdjustVariableBounds(const Vector& new_x_L,
203                                       const Vector& new_x_U,
204                                       const Vector& new_d_L,
205                                       const Vector& new_d_U);
206 
207     /** @name Counters for the number of function evaluations. */
208     //@{
f_evals() const209     virtual Index f_evals() const
210     {
211       return f_evals_;
212     }
grad_f_evals() const213     virtual Index grad_f_evals() const
214     {
215       return grad_f_evals_;
216     }
c_evals() const217     virtual Index c_evals() const
218     {
219       return c_evals_;
220     }
jac_c_evals() const221     virtual Index jac_c_evals() const
222     {
223       return jac_c_evals_;
224     }
d_evals() const225     virtual Index d_evals() const
226     {
227       return d_evals_;
228     }
jac_d_evals() const229     virtual Index jac_d_evals() const
230     {
231       return jac_d_evals_;
232     }
h_evals() const233     virtual Index h_evals() const
234     {
235       return h_evals_;
236     }
237     //@}
238 
239     /** Solution Routines - overloaded from IpoptNLP*/
240     //@{
241     void FinalizeSolution(SolverReturn status,
242                           const Vector& x, const Vector& z_L, const Vector& z_U,
243                           const Vector& c, const Vector& d,
244                           const Vector& y_c, const Vector& y_d,
245                           Number obj_value,
246                           const IpoptData* ip_data,
247                           IpoptCalculatedQuantities* ip_cq);
248     bool IntermediateCallBack(AlgorithmMode mode,
249                               Index iter, Number obj_value,
250                               Number inf_pr, Number inf_du,
251                               Number mu, Number d_norm,
252                               Number regularization_size,
253                               Number alpha_du, Number alpha_pr,
254                               Index ls_trials,
255                               SmartPtr<const IpoptData> ip_data,
256                               SmartPtr<IpoptCalculatedQuantities> ip_cq);
257     //@}
258 
259     /** @name Methods for IpoptType */
260     //@{
261     /** Called by IpoptType to register the options */
262     static void RegisterOptions(SmartPtr<RegisteredOptions> roptions);
263     //@}
264 
265     /** Accessor method to the underlying NLP */
nlp()266     SmartPtr<NLP> nlp()
267     {
268       return nlp_;
269     }
270 
271     /**@name Methods related to function evaluation timing. */
272     //@{
273 
274     /** Reset the timing statistics */
275     void ResetTimes();
276 
277     void PrintTimingStatistics(Journalist& jnlst,
278                                EJournalLevel level,
279                                EJournalCategory category) const;
280 
f_eval_time() const281     const TimedTask& f_eval_time() const
282     {
283       return f_eval_time_;
284     }
grad_f_eval_time() const285     const TimedTask& grad_f_eval_time() const
286     {
287       return grad_f_eval_time_;
288     }
c_eval_time() const289     const TimedTask& c_eval_time() const
290     {
291       return c_eval_time_;
292     }
jac_c_eval_time() const293     const TimedTask& jac_c_eval_time() const
294     {
295       return jac_c_eval_time_;
296     }
d_eval_time() const297     const TimedTask& d_eval_time() const
298     {
299       return d_eval_time_;
300     }
jac_d_eval_time() const301     const TimedTask& jac_d_eval_time() const
302     {
303       return jac_d_eval_time_;
304     }
h_eval_time() const305     const TimedTask& h_eval_time() const
306     {
307       return h_eval_time_;
308     }
309 
310     Number TotalFunctionEvaluationCpuTime() const;
311     Number TotalFunctionEvaluationSysTime() const;
312     Number TotalFunctionEvaluationWallclockTime() const;
313     //@}
314 
315   private:
316     /** journalist */
317     SmartPtr<const Journalist> jnlst_;
318 
319     /** Pointer to the NLP */
320     SmartPtr<NLP> nlp_;
321 
322     /** Necessary Vector/Matrix spaces */
323     //@{
324     SmartPtr<const VectorSpace> x_space_;
325     SmartPtr<const VectorSpace> c_space_;
326     SmartPtr<const VectorSpace> d_space_;
327     SmartPtr<const VectorSpace> x_l_space_;
328     SmartPtr<const MatrixSpace> px_l_space_;
329     SmartPtr<const VectorSpace> x_u_space_;
330     SmartPtr<const MatrixSpace> px_u_space_;
331     SmartPtr<const VectorSpace> d_l_space_;
332     SmartPtr<const MatrixSpace> pd_l_space_;
333     SmartPtr<const VectorSpace> d_u_space_;
334     SmartPtr<const MatrixSpace> pd_u_space_;
335     SmartPtr<const MatrixSpace> jac_c_space_;
336     SmartPtr<const MatrixSpace> jac_d_space_;
337     SmartPtr<const SymMatrixSpace> h_space_;
338 
339     SmartPtr<const MatrixSpace> scaled_jac_c_space_;
340     SmartPtr<const MatrixSpace> scaled_jac_d_space_;
341     SmartPtr<const SymMatrixSpace> scaled_h_space_;
342     //@}
343     /**@name Storage for Model Quantities */
344     //@{
345     /** Objective function */
346     CachedResults<Number> f_cache_;
347 
348     /** Gradient of the objective function */
349     CachedResults<SmartPtr<const Vector> > grad_f_cache_;
350 
351     /** Equality constraint residuals */
352     CachedResults<SmartPtr<const Vector> > c_cache_;
353 
354     /** Jacobian Matrix for equality constraints
355      *  (current iteration) */
356     CachedResults<SmartPtr<const Matrix> > jac_c_cache_;
357 
358     /** Inequality constraint residual (reformulated
359      *  as equalities with slacks */
360     CachedResults<SmartPtr<const Vector> > d_cache_;
361 
362     /** Jacobian Matrix for inequality constraints
363      *  (current iteration) */
364     CachedResults<SmartPtr<const Matrix> > jac_d_cache_;
365 
366     /** Hessian of the lagrangian
367      *  (current iteration) */
368     CachedResults<SmartPtr<const SymMatrix> > h_cache_;
369 
370     /** Unscaled version of x vector */
371     CachedResults<SmartPtr<const Vector> > unscaled_x_cache_;
372 
373     /** Lower bounds on x */
374     SmartPtr<const Vector> x_L_;
375 
376     /** Permutation matrix (x_L_ -> x) */
377     SmartPtr<const Matrix> Px_L_;
378 
379     /** Upper bounds on x */
380     SmartPtr<const Vector> x_U_;
381 
382     /** Permutation matrix (x_U_ -> x */
383     SmartPtr<const Matrix> Px_U_;
384 
385     /** Lower bounds on d */
386     SmartPtr<const Vector> d_L_;
387 
388     /** Permutation matrix (d_L_ -> d) */
389     SmartPtr<const Matrix> Pd_L_;
390 
391     /** Upper bounds on d */
392     SmartPtr<const Vector> d_U_;
393 
394     /** Permutation matrix (d_U_ -> d */
395     SmartPtr<const Matrix> Pd_U_;
396 
397     /** Original unmodified lower bounds on x */
398     SmartPtr<const Vector> orig_x_L_;
399 
400     /** Original unmodified upper bounds on x */
401     SmartPtr<const Vector> orig_x_U_;
402     //@}
403 
404     /**@name Default Compiler Generated Methods
405      * (Hidden to avoid implicit creation/calling).
406      * These methods are not implemented and
407      * we do not want the compiler to implement
408      * them for us, so we declare them private
409      * and do not define them. This ensures that
410      * they will not be implicitly created/called. */
411     //@{
412     /** Default Constructor */
413     OrigIpoptNLP();
414 
415     /** Copy Constructor */
416     OrigIpoptNLP(const OrigIpoptNLP&);
417 
418     /** Overloaded Equals Operator */
419     void operator=(const OrigIpoptNLP&);
420     //@}
421 
422     /** @name auxilliary functions */
423     //@{
424     /** relax the bounds by a relative move of relax_bound_factor.
425      *  Here, relax_bound_factor should be negative (or zero) for
426      *  lower bounds, and positive (or zero) for upper bounds.
427      */
428     void relax_bounds(Number bound_relax_factor, Vector& bounds);
429     /** Method for getting the unscaled version of the x vector */
430     SmartPtr<const Vector> get_unscaled_x(const Vector& x);
431     //@}
432 
433     /** @name Algorithmic parameters */
434     //@{
435     /** relaxation factor for the bounds */
436     Number bound_relax_factor_;
437     /** Flag indicating whether the primal variables should be
438      *  projected back into original bounds are optimization. */
439     bool honor_original_bounds_;
440     /** Flag indicating whether the TNLP with identical structure has
441      *  already been solved before. */
442     bool warm_start_same_structure_;
443     /** Flag indicating what Hessian information is to be used. */
444     HessianApproximationType hessian_approximation_;
445     /** Flag indicating in which space Hessian is to be approximated. */
446     HessianApproximationSpace hessian_approximation_space_;
447     /** Flag indicating whether it is desired to check if there are
448      *  Nan or Inf entries in first and second derivative matrices. */
449     bool check_derivatives_for_naninf_;
450     /** Flag indicating if we need to ask for equality constraint
451      *  Jacobians only once */
452     bool jac_c_constant_;
453     /** Flag indicating if we need to ask for inequality constraint
454      *  Jacobians only once */
455     bool jac_d_constant_;
456     /** Flag indicating if we need to ask for Hessian only once */
457     bool hessian_constant_;
458     //@}
459 
460     /** @name Counters for the function evaluations */
461     //@{
462     Index f_evals_;
463     Index grad_f_evals_;
464     Index c_evals_;
465     Index jac_c_evals_;
466     Index d_evals_;
467     Index jac_d_evals_;
468     Index h_evals_;
469     //@}
470 
471     /** Flag indicating if initialization method has been called */
472     bool initialized_;
473 
474     /**@name Timing statistics for the function evaluations. */
475     //@{
476     TimedTask f_eval_time_;
477     TimedTask grad_f_eval_time_;
478     TimedTask c_eval_time_;
479     TimedTask jac_c_eval_time_;
480     TimedTask d_eval_time_;
481     TimedTask jac_d_eval_time_;
482     TimedTask h_eval_time_;
483     //@}
484   };
485 
486 } // namespace Ipopt
487 
488 #endif
489