CppAD: A C++ Algorithmic Differentiation Package  20130102
cppad_ipopt_nlp.hpp
Go to the documentation of this file.
00001 /* $Id$ */
00002 # ifndef CPPAD_CPPAD_IPOPT_NLP_INCLUDED
00003 # define CPPAD_CPPAD_IPOPT_NLP_INCLUDED
00004 /* --------------------------------------------------------------------------
00005 CppAD: C++ Algorithmic Differentiation: Copyright (C) 2003-12 Bradley M. Bell
00006 
00007 CppAD is distributed under multiple licenses. This distribution is under
00008 the terms of the 
00009                     Eclipse Public License Version 1.0.
00010 
00011 A copy of this license is included in the COPYING file of this distribution.
00012 Please visit http://www.coin-or.org/CppAD/ for information on other licenses.
00013 -------------------------------------------------------------------------- */
00014 /*
00015 $begin cppad_ipopt_nlp$$
00016 $dollar @$$
00017 $spell
00018      libipopt
00019      namespace
00020      dir
00021      cppad
00022      bool
00023      doesn't
00024      nan
00025      inf
00026      naninf
00027      std
00028      maxiter
00029      infeasibility
00030      obj
00031      const
00032      optimizer
00033      cppad_ipopt_nlp.hpp
00034      fg_info.eval
00035      retape
00036      CppAD
00037      
00038 $$
00039 $section Nonlinear Programming Using the CppAD Interface to Ipopt$$
00040 
00041 $index nonlinear, programming CppAD$$
00042 $index deprecated, ipopt interface$$
00043 $index programming, nonlinear$$
00044 $index CppAD, nonlinear programming$$
00045 $index Ipopt, AD$$
00046 $index AD, Ipopt$$
00047 
00048 $head Deprecated$$
00049 This interface to Ipopt is deprecated, use $cref ipopt_solve$$ instead.
00050 
00051 $head Syntax$$
00052 $codei%# include "cppad_ipopt_nlp.hpp"
00053 %$$
00054 $codei%cppad_ipopt_solution %solution%;
00055 %$$
00056 $codei%cppad_ipopt_nlp %cppad_nlp%(
00057      %n%, %m%, %x_i%, %x_l%, %x_u%, %g_l%, %g_u%, &%fg_info%, &%solution%
00058 )%$$
00059 $codei%
00060 export LD_LIBRARY_PATH=@LD_LIBRARY_PATH:%ipopt_library_paths%$$
00061 
00062 $head Purpose$$
00063 The class $code cppad_ipopt_nlp$$ is used to solve nonlinear programming
00064 problems of the form
00065 $latex \[
00066 \begin{array}{rll}
00067 {\rm minimize}      & f(x) 
00068 \\
00069 {\rm subject \; to} & g^l \leq g(x) \leq g^u
00070 \\
00071                     & x^l  \leq x   \leq x^u
00072 \end{array}
00073 \] $$
00074 This is done using 
00075 $href%
00076      http://www.coin-or.org/projects/Ipopt.xml%
00077      Ipopt
00078 %$$
00079 optimizer and 
00080 $href%
00081      http://www.coin-or.org/CppAD/%
00082      CppAD
00083 %$$
00084 Algorithmic Differentiation package.
00085 
00086 $head cppad_ipopt namespace$$
00087 All of the declarations for these routines
00088 are in the $code cppad_ipopt$$ namespace 
00089 (not the $code CppAD$$ namespace).
00090 For example; $cref/SizeVector/cppad_ipopt_nlp/SizeVector/$$ below
00091 actually denotes the type $code cppad_ipopt::SizeVector$$.
00092 
00093 $head ipopt_library_paths$$
00094 If you are linking to a shared version of the Ipopt library,
00095 you may have to add some paths the $code LD_LIBRARY_PATH$$ 
00096 shell variable using the $code export$$ command in the syntax above.
00097 For example, if the file the ipopt library is
00098 $codei%
00099      %ipopt_prefix%/lib64/libipopt.a
00100 %$$
00101 you will need to add the corresponding directory; e.g.,
00102 $code%
00103      export LD_LIBRARY_PATH="%ipopt_prefix%/lib64%:@LD_LIBRARY_PATH"
00104 %$$
00105 see $cref ipopt_prefix$$.
00106 
00107 $head fg(x)$$
00108 The function $latex fg : \B{R}^n \rightarrow \B{R}^{m+1}$$ is defined by
00109 $latex \[
00110 \begin{array}{rcl}
00111      fg_0 (x)     & = & f(x)         \\
00112      fg_1 (x)     & = & g_0 (x)      \\
00113                   & \vdots &         \\
00114      fg_m (x)     & = & g_{m-1} (x)
00115      \end{array}
00116 \] $$
00117 
00118 $subhead Index Vector$$
00119 We define an $icode index vector$$ as a vector of non-negative integers
00120 for which none of the values are equal; i.e.,
00121 it is both a vector and a set.
00122 If $latex I$$ is an index vector $latex |I|$$ is used to denote the
00123 number of elements in $latex I$$ and $latex \| I \|$$ is used
00124 to denote the value of the maximum element in $latex I$$.
00125 
00126 $subhead Projection$$
00127 Given an index vector $latex J$$ and a positive integer $latex n$$
00128 where $latex n > \| J \|$$, we use $latex J \otimes n $$ for
00129 the mapping $latex ( J \otimes n ) : \B{R}^n \rightarrow \B{R}^{|J|}$$ defined by
00130 $latex \[
00131      [ J \otimes n ] (x)_j = x_{J(j)}
00132 \] $$
00133 for $latex j = 0 , \ldots |J| - 1$$.
00134 
00135 $subhead Injection$$
00136 Given an index vector $latex I$$ and a positive integer $latex m$$
00137 where $latex m > \| I \|$$, we use $latex m \otimes I$$ for
00138 the mapping $latex ( m \otimes I ): \B{R}^{|I|} \rightarrow \B{R}^m$$ defined by
00139 $latex \[
00140 [ m \otimes I ] (y)_i = \left\{ \begin{array}{ll}
00141 y_k & {\rm if} \; i = I(k) \; {\rm for \; some} \; 
00142      k \in \{ 0 , \cdots, |I|-1 \} 
00143 \\
00144 0   & {\rm otherwise}
00145 \end{array} \right.
00146 \] $$
00147 
00148 $subhead Representation$$
00149 In many applications, each of the component functions of $latex fg(x)$$
00150 only depend on a few of the components of $latex x$$.
00151 In this case, expressing $latex fg(x)$$ in terms of simpler functions
00152 with fewer arguments can greatly reduce the amount of work required
00153 to compute its derivatives.
00154 $pre
00155 
00156 $$
00157 We use the functions
00158 $latex r_k : \B{R}^{q(k)} \rightarrow \B{R}^{p(k)}$$ 
00159 for $latex k = 0 , \ldots , K$$ to express our
00160 representation of $latex fg(x)$$ in terms of simpler functions
00161 as follows
00162 $latex \[
00163 fg(x) = \sum_{k=0}^{K-1} \; \sum_{\ell=0}^{L(k) - 1} 
00164 [ (m+1) \otimes I_{k,\ell} ] \; \circ
00165       \; r_k \; \circ \; [ J_{k,\ell} \otimes n ] \; (x)
00166 \] $$
00167 where $latex \circ$$ represents function composition,
00168 for $latex k = 0 , \ldots , K - 1$$, and $latex \ell = 0 , \ldots , L(k)$$,
00169 $latex I_{k,\ell}$$ and  $latex J_{k,\ell}$$ are index vectors with
00170 $latex | J_{k,\ell} | = q(k)$$, 
00171 $latex \| J_{k,\ell} \| < n$$,
00172 $latex | I_{k,\ell} | = p(k)$$, and
00173 $latex \| I_{k,\ell} \| \leq m$$. 
00174 
00175 $head Simple Representation$$
00176 In the simple representation,
00177 $latex r_0 (x) = fg(x)$$,
00178 $latex K = 1$$,
00179 $latex q(0) = n$$,
00180 $latex p(0) = m+1$$,
00181 $latex L(0) = 1$$,
00182 $latex I_{0,0} = (0 , \ldots , m)$$,
00183 and $latex J_{0,0} = (0 , \ldots , n-1)$$.
00184 
00185 $head SizeVector$$
00186 The type $codei SizeVector$$ is defined by the 
00187 $codei cppad_ipopt_nlp.hpp$$ include file to be a 
00188 $cref SimpleVector$$ class with elements of type
00189 $code size_t$$.
00190 
00191 $head NumberVector$$
00192 The type $codei NumberVector$$ is defined by the 
00193 $codei cppad_ipopt_nlp.hpp$$ include file to be a 
00194 $cref SimpleVector$$ class with elements of type
00195 $code Ipopt::Number$$.
00196 
00197 $head ADNumber$$
00198 The type $codei ADNumber$$ is defined by the 
00199 $codei cppad_ipopt_nlp.hpp$$ include file to be a 
00200 an AD type that can be used to compute derivatives.
00201 
00202 $head ADVector$$
00203 The type $codei ADVector$$ is defined by the 
00204 $codei cppad_ipopt_nlp.hpp$$ include file to be a 
00205 $cref SimpleVector$$ class with elements of type
00206 $code ADNumber$$. 
00207 
00208 $head n$$
00209 The argument $icode n$$ has prototype
00210 $codei%
00211      size_t %n%
00212 %$$
00213 It specifies the dimension of the argument space; 
00214 i.e., $latex x \in \B{R}^n$$.
00215 
00216 $head m$$
00217 The argument $icode m$$ has prototype
00218 $codei%
00219      size_t %m%
00220 %$$
00221 It specifies the dimension of the range space for $latex g$$; 
00222 i.e., $latex g : \B{R}^n \rightarrow \B{R}^m$$.
00223 
00224 $head x_i$$
00225 The argument $icode x_i$$ has prototype
00226 $codei%
00227      const NumberVector& %x_i%
00228 %$$
00229 and its size is equal to $latex n$$.
00230 It specifies the initial point where Ipopt starts the optimization process.
00231 
00232 $head x_l$$
00233 The argument $icode x_l$$ has prototype
00234 $codei%
00235      const NumberVector& %x_l%
00236 %$$
00237 and its size is equal to $latex n$$.
00238 It specifies the lower limits for the argument in the optimization problem;
00239 i.e., $latex x^l$$.
00240 
00241 $head x_u$$
00242 The argument $icode x_u$$ has prototype
00243 $codei%
00244      const NumberVector& %x_u%
00245 %$$
00246 and its size is equal to $latex n$$.
00247 It specifies the upper limits for the argument in the optimization problem;
00248 i.e., $latex x^u$$.
00249 
00250 $head g_l$$
00251 The argument $icode g_l$$ has prototype
00252 $codei%
00253      const NumberVector& %g_l%
00254 %$$
00255 and its size is equal to $latex m$$.
00256 It specifies the lower limits for the constraints in the optimization problem;
00257 i.e., $latex g^l$$.
00258 
00259 $head g_u$$
00260 The argument $icode g_u$$ has prototype
00261 $codei%
00262      const NumberVector& %g_u%
00263 %$$
00264 and its size is equal to $latex n$$.
00265 It specifies the upper limits for the constraints in the optimization problem;
00266 i.e., $latex g^u$$.
00267 
00268 $head fg_info$$
00269 The argument $icode fg_info$$ has prototype
00270 $codei%
00271      %FG_info fg_info%
00272 %$$
00273 where the class $icode FG_info$$ is derived from the 
00274 base class $code cppad_ipopt_fg_info$$.
00275 Certain virtual member functions of $icode fg_info$$ are used to 
00276 compute the value of $latex fg(x)$$.
00277 The specifications for these member functions are given below:
00278 
00279 $subhead fg_info.number_functions$$
00280 This member function has prototype
00281 $codei%
00282      virtual size_t cppad_ipopt_fg_info::number_functions(void)
00283 %$$
00284 If $icode K$$ has type $code size_t$$, the syntax
00285 $codei%
00286      %K% = %fg_info%.number_functions()
00287 %$$
00288 sets $icode K$$ to the number of functions used in the
00289 representation of $latex fg(x)$$; i.e., $latex K$$ in
00290 the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00291 $pre
00292 
00293 $$
00294 The $code cppad_ipopt_fg_info$$ implementation of this function
00295 corresponds to the simple representation mentioned above; i.e.
00296 $icode%K% = 1%$$.
00297 
00298 $subhead fg_info.eval_r$$
00299 This member function has the prototype
00300 $codei%
00301 virtual ADVector cppad_ipopt_fg_info::eval_r(size_t %k%, const ADVector& %u%) = 0;
00302 %$$
00303 Thus it is a pure virtual function and must be defined in the 
00304 derived class $icode FG_info$$.
00305 $pre
00306 
00307 $$
00308 This function computes the value of $latex r_k (u)$$
00309 used in the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$
00310 for $latex fg(x)$$.
00311 If $icode k$$ in $latex \{0 , \ldots , K-1 \}$$ has type $code size_t$$,
00312 $icode u$$ is an $code ADVector$$ of size $icode q(k)$$ 
00313 and $icode r$$ is an $code ADVector$$ of size $icode p(k)$$
00314 the syntax
00315 $codei%
00316      %r% = %fg_info%.eval_r(%k%, %u%)
00317 %$$
00318 set $icode r$$ to the vector $latex r_k (u)$$.
00319 
00320 $subhead fg_info.retape$$
00321 This member function has the prototype
00322 $codei%
00323      virtual bool cppad_ipopt_fg_info::retape(size_t %k%)
00324 %$$
00325 If $icode k$$ in $latex \{0 , \ldots , K-1 \}$$ has type $code size_t$$,
00326 and $icode retape$$ has type $code bool$$,
00327 the syntax
00328 $codei%
00329         %retape% = %fg_info%.retape(%k%)
00330 %$$
00331 sets $icode retape$$ to true or false.
00332 If $icode retape$$ is true, 
00333 $code cppad_ipopt_nlp$$ will retape the operation sequence 
00334 corresponding to $latex r_k (u)$$ for
00335 every value of $icode u$$. 
00336 An $code cppad_ipopt_nlp$$ object
00337 should use much less memory and run faster if $icode retape$$ is false.
00338 You can test both the true and false cases to make sure 
00339 the operation sequence does not depend on $icode u$$.
00340 $pre
00341 
00342 $$
00343 The $code cppad_ipopt_fg_info$$ implementation of this function
00344 sets $icode retape$$ to true 
00345 (while slower it is also safer to always retape).
00346 
00347 $subhead fg_info.domain_size$$
00348 This member function has prototype
00349 $codei%
00350      virtual size_t cppad_ipopt_fg_info::domain_size(size_t %k%)
00351 %$$
00352 If $icode k$$ in $latex \{0 , \ldots , K-1 \}$$ has type $code size_t$$,
00353 and $icode q$$ has type $code size_t$$, the syntax
00354 $codei%
00355      %q% = %fg_info%.domain_size(%k%)
00356 %$$
00357 sets $icode q$$ to the dimension of the domain space for $latex r_k (u)$$;
00358 i.e., $latex q(k)$$ in
00359 the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00360 
00361 $pre
00362 
00363 $$
00364 The $code cppad_ipopt_h_base$$ implementation of this function
00365 corresponds to the simple representation mentioned above; i.e.,
00366 $latex q = n$$.
00367 
00368 $subhead fg_info.range_size$$
00369 This member function has prototype
00370 $codei%
00371      virtual size_t cppad_ipopt_fg_info::range_size(size_t %k%)
00372 %$$
00373 If $icode k$$ in $latex \{0 , \ldots , K-1 \}$$ has type $code size_t$$,
00374 and $icode p$$ has type $code size_t$$, the syntax
00375 $codei%
00376      %p% = %fg_info%.range_size(%k%)
00377 %$$
00378 sets $icode p$$ to the dimension of the range space for $latex r_k (u)$$;
00379 i.e., $latex p(k)$$ in
00380 the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00381 $pre
00382 
00383 $$
00384 The $code cppad_ipopt_h_base$$ implementation of this function
00385 corresponds to the simple representation mentioned above; i.e.,
00386 $latex p = m+1$$.
00387 
00388 $subhead fg_info.number_terms$$
00389 This member function has prototype
00390 $codei%
00391      virtual size_t cppad_ipopt_fg_info::number_terms(size_t %k%)
00392 %$$
00393 If $icode k$$ in $latex \{0 , \ldots , K-1 \}$$ has type $code size_t$$,
00394 and $icode L$$ has type $code size_t$$, the syntax
00395 $codei%
00396      %L% = %fg_info%.number_terms(%k%)
00397 %$$
00398 sets $icode L$$ to the number of terms in representation
00399 for this value of $icode k$$;
00400 i.e., $latex L(k)$$ in
00401 the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00402 $pre
00403 
00404 $$
00405 The $code cppad_ipopt_h_base$$ implementation of this function
00406 corresponds to the simple representation mentioned above; i.e.,
00407 $latex L = 1$$.
00408 
00409 $subhead fg_info.index$$
00410 This member function has prototype
00411 $codei%
00412      virtual void cppad_ipopt_fg_info::index(
00413           size_t %k%, size_t %ell%, SizeVector& %I%, SizeVector& %J%
00414      )
00415 %$$ 
00416 The argument 
00417 $icode% 
00418      k
00419 %$$
00420 has type $codei size_t$$
00421 and is a value between zero and $latex K-1$$ inclusive.
00422 The argument 
00423 $icode% 
00424      ell
00425 %$$
00426 has type $codei size_t$$
00427 and is a value between zero and $latex L(k)-1$$ inclusive.
00428 The argument 
00429 $icode%
00430      I
00431 %$$ is a $cref SimpleVector$$ with elements
00432 of type $code size_t$$ and size greater than or equal to $latex p(k)$$.
00433 The input value of the elements of $icode I$$ does not matter.
00434 The output value of
00435 the first $latex p(k)$$ elements of $icode I$$ 
00436 must be the corresponding elements of $latex I_{k,ell}$$ 
00437 in the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00438 The argument 
00439 $icode%
00440      J
00441 %$$ is a $cref SimpleVector$$ with elements
00442 of type $code size_t$$ and size greater than or equal to $latex q(k)$$.
00443 The input value of the elements of $icode J$$ does not matter.
00444 The output value of 
00445 the first $latex q(k)$$ elements of $icode J$$ 
00446 must be the corresponding elements of $latex J_{k,ell}$$ 
00447 in the $cref/representation/cppad_ipopt_nlp/fg(x)/Representation/$$ above.
00448 $pre
00449 
00450 $$
00451 The $code cppad_ipopt_h_base$$ implementation of this function
00452 corresponds to the simple representation mentioned above; i.e.,
00453 for $latex i = 0 , \ldots , m$$,
00454 $icode%I%[%i%] = %i%$$,
00455 and  for $latex j = 0 , \ldots , n-1$$,
00456 $icode%J%[%j%] = %j%$$.
00457 
00458 $head solution$$
00459 After the optimization process is completed, $icode solution$$ contains
00460 the following information:
00461 
00462 $subhead status$$
00463 The $icode status$$ field of $icode solution$$ has prototype
00464 $codei%
00465      cppad_ipopt_solution::solution_status %solution%.status
00466 %$$
00467 It is the final Ipopt status for the optimizer. 
00468 Here is a list of the possible values for the status:
00469 
00470 $table
00471 $icode status$$ $cnext Meaning
00472 $rnext
00473 not_defined $cnext
00474 The optimizer did not return a final status to this $code cppad_ipopt_nlp$$
00475 object.
00476 $rnext
00477 unknown $cnext
00478 The status returned by the optimizer is not defined in the Ipopt
00479 documentation for $code finalize_solution$$.
00480 $rnext
00481 success $cnext
00482 Algorithm terminated successfully at a point satisfying the convergence 
00483 tolerances (see Ipopt options).
00484 $rnext
00485 maxiter_exceeded $cnext
00486 The maximum number of iterations was exceeded (see Ipopt options).
00487 $rnext
00488 stop_at_tiny_step $cnext
00489 Algorithm terminated because progress was very slow.
00490 $rnext
00491 stop_at_acceptable_point $cnext
00492 Algorithm stopped at a point that was converged, 
00493 not to the 'desired' tolerances, but to 'acceptable' tolerances 
00494 (see Ipopt options).
00495 $rnext
00496 local_infeasibility $cnext
00497 Algorithm converged to a non-feasible point
00498 (problem may have no solution).
00499 $rnext
00500 user_requested_stop $cnext
00501 This return value should not happen.
00502 $rnext
00503 diverging_iterates $cnext
00504 It the iterates are diverging.
00505 $rnext
00506 restoration_failure $cnext
00507 Restoration phase failed, algorithm doesn't know how to proceed.
00508 $rnext
00509 error_in_step_computation $cnext
00510 An unrecoverable error occurred while Ipopt tried to 
00511 compute the search direction.
00512 $rnext
00513 invalid_number_detected $cnext
00514 Algorithm received an invalid number (such as $code nan$$ or $code inf$$) 
00515 from the users function $icode%fg_info%.eval%$$ or from the CppAD evaluations
00516 of its derivatives
00517 (see the Ipopt option $code check_derivatives_for_naninf$$).
00518 $rnext
00519 internal_error $cnext
00520 An unknown Ipopt internal error occurred.
00521 Contact the Ipopt authors through the mailing list.
00522 $tend
00523 
00524 $subhead x$$
00525 The $code x$$ field of $icode solution$$ has prototype
00526 $codei%
00527      NumberVector %solution%.x
00528 %$$
00529 and its size is equal to $latex n$$.
00530 It is the final $latex x$$ value for the optimizer.
00531 
00532 $subhead z_l$$
00533 The $code z_l$$ field of $icode solution$$ has prototype
00534 $codei%
00535      NumberVector %solution%.z_l
00536 %$$
00537 and its size is equal to $latex n$$.
00538 It is the final Lagrange multipliers for the 
00539 lower bounds on $latex x$$.
00540 
00541 $subhead z_u$$
00542 The $code z_u$$ field of $icode solution$$ has prototype
00543 $codei%
00544      NumberVector %solution%.z_u
00545 %$$
00546 and its size is equal to $latex n$$.
00547 It is the final Lagrange multipliers for the 
00548 upper bounds on $latex x$$.
00549 
00550 $subhead g$$
00551 The $code g$$ field of $icode solution$$ has prototype
00552 $codei%
00553      NumberVector %solution%.g
00554 %$$
00555 and its size is equal to $latex m$$.
00556 It is the final value for the constraint function $latex g(x)$$.
00557 
00558 $subhead lambda$$
00559 The $code lambda$$ field of $icode solution$$ has prototype
00560 $codei%
00561      NumberVector %solution%.lambda
00562 %$$
00563 and its size is equal to $latex m$$.
00564 It is the final value for the 
00565 Lagrange multipliers corresponding to the constraint function.
00566 
00567 $subhead obj_value$$
00568 The $code obj_value$$ field of $icode solution$$ has prototype
00569 $codei%
00570      Number %solution%.obj_value
00571 %$$
00572 It is the final value of the objective function $latex f(x)$$.
00573 
00574 
00575 $children%
00576      cppad_ipopt/example/get_started.cpp%
00577      cppad_ipopt/example/ode1.omh%
00578      cppad_ipopt/speed/ode_speed.cpp
00579 %$$
00580 
00581 $head Example$$
00582 The file 
00583 $cref ipopt_nlp_get_started.cpp$$ is an example and test of 
00584 $code cppad_ipopt_nlp$$  that uses the 
00585 $cref/simple representation/cppad_ipopt_nlp/Simple Representation/$$.
00586 It returns true if it succeeds and false otherwise.
00587 The section $cref ipopt_nlp_ode$$ discusses an example that
00588 uses a more complex representation.
00589 
00590 $head Wish List$$
00591 This is a list of possible future improvements to 
00592 $code cppad_ipopt_nlp$$ that would require changed to the user interface:
00593 $list number$$
00594 The routine $codei%fg_info.eval_r(%k%, %u%)%$$ should also support 
00595 $codei NumberVector$$ for the type of the argument $code u$$
00596 (this would certainly be more efficient when 
00597 $codei%fg_info.retape(%k%)%$$ is true and $latex L(k) > 1$$).
00598 It could be an option for the user to provide this as well as
00599 the necessary $code ADVector$$ definition.
00600 $lnext
00601 There should a $cref Discrete$$ routine that the user can call
00602 to determine the value of $latex \ell$$ during the evaluation of
00603 $codei%fg_info.eval_r(%k%, %u%)%$$.
00604 This way data, which does not affect the derivative values,
00605 can be included in the function recording and evaluation.
00606 $lend
00607 
00608 $end
00609 -----------------------------------------------------------------------------
00610 */
00611 # include <cppad/cppad.hpp>
00612 # include <coin/IpIpoptApplication.hpp>
00613 # include <coin/IpTNLP.hpp>
00614 
00615 /*!
00616 \file cppad_ipopt_nlp.hpp
00617 \brief CppAD interface to Ipopt
00618 
00619 \ingroup cppad_ipopt_nlp_cpp
00620 */
00621 
00622 // ---------------------------------------------------------------------------
00623 namespace cppad_ipopt {
00624 // ---------------------------------------------------------------------------
00625 
00626 /// A scalar value used to record operation sequence.
00627 typedef CppAD::AD<Ipopt::Number>       ADNumber;
00628 /// A simple vector of values used to record operation sequence
00629 typedef CppAD::vector<ADNumber>        ADVector;
00630 /// A simple vector of size_t values.
00631 typedef CppAD::vector<size_t>          SizeVector;
00632 /// A simple vector of values used by Ipopt
00633 typedef CppAD::vector<Ipopt::Number>   NumberVector;
00634 
00635 /*!
00636 Abstract base class user derives from to define the funcitons in the problem.
00637 */
00638 class cppad_ipopt_fg_info
00639 {
00640      /// allow cppad_ipopt_nlp class complete access to this class
00641      friend class cppad_ipopt_nlp;
00642 private:
00643      /// domain space dimension for the functions f(x), g(x)
00644      size_t n_;
00645      /// range space dimension for the function g(x)
00646      size_t m_;
00647      /// the cppad_ipopt_nlp constructor uses this method to set n_
00648      void set_n(size_t n)
00649      {    n_ = n; }
00650      /// the cppad_ipopt_nlp constructor uses this method to set m_
00651      void set_m(size_t m)
00652      {    m_ = m; }
00653 
00654 public:
00655      /// destructor virtual so user derived class destructor gets called
00656      virtual ~cppad_ipopt_fg_info(void)
00657      { }
00658      /// number_functions; i.e. K (simple representation uses 1)
00659      virtual size_t number_functions(void)
00660      {    return 1; }
00661      /// function that evaluates the users representation for f(x) and
00662      /// and g(x) is pure virtual so user must define it in derived class
00663      virtual ADVector eval_r(size_t k, const ADVector& u) = 0;
00664      /// should the function r_k (u) be retaped when ever the arguemnt
00665      /// u changes (default is true which is safe but slow)
00666      virtual bool retape(size_t k)
00667      {    return true; }
00668      /// domain_size q[k] for r_k (u) (simple representation uses n)
00669      virtual size_t domain_size(size_t k)
00670      {    return n_; }
00671      /// range_size p[k] for r_k (u) (simple representation uses m+1)
00672      virtual size_t range_size(size_t k)
00673      {    return m_ + 1; }
00674      /// number_terms that use r_k (u) (simple represenation uses 1)
00675      virtual size_t number_terms(size_t k)
00676      {    return 1; }
00677      /// return the index vectors I_{k,ell} and J_{k,ell}
00678      /// (simple representation uses I[i] = i and J[j] = j)
00679      virtual void index(size_t k, size_t ell, SizeVector& I, SizeVector& J)
00680      {    assert( I.size() >= m_ + 1 );
00681           assert( J.size() >= n_ );
00682           for(size_t i = 0; i <= m_; i++)
00683                I[i] = i;
00684           for(size_t j = 0; j < n_; j++)
00685                J[j] = j;
00686      }
00687 };
00688 
00689 /*!
00690 Class that contains information about the problem solution
00691 
00692 \section Nonlinear_Programming_Problem Nonlinear Programming Problem
00693 We are give smooth functions 
00694 \f$ f : {\bf R}^n \rightarrow {\bf R} \f$
00695 and
00696 \f$ g : {\bf R}^n \rightarrow {\bf R}^m \f$
00697 and wish to solve the problem
00698 \f[
00699 \begin{array}{rcl}
00700 {\rm minimize} & f(x) & {\rm w.r.t.} \; x \in {\bf R}^n
00701 \\
00702 {\rm subject \; to} & g^l \leq g(x) \leq g^u 
00703 \\
00704 & x^l \leq x \leq x^u
00705 \end{array}
00706 \f]
00707 
00708 
00709 \section Users_Representation Users Representation
00710 The functions 
00711 \f$ f : {\bf R}^n \rightarrow {\bf R} \f$ and
00712 \f$ g : {\bf R}^n \rightarrow {\bf R}^m \f$ are defined by
00713 \f[
00714 \left( \begin{array}{c} f(x) \\ g(x) \end{array} \right)
00715 =
00716 \sum_{k=0}^{K-1} \; \sum_{\ell=0}^{L(k) - 1} 
00717 [ (m+1) \otimes I_{k,\ell} ] \; \circ
00718       \; r_k \; \circ \; [ J_{k,\ell} \otimes n ] \; (x)
00719 \f]
00720 where for \f$ k = 0 , \ldots , K-1\f$,
00721 \f$ r_k : {\bf R}^{q(k)} \rightarrow {\bf R}^{p(k)} \f$.
00722 
00723 \section Evaluation_Methods Evaluation Methods
00724 The set of evaluation methods for this class is
00725 \verbatim
00726      { eval_f, eval_grad_f, eval_g, eval_jac_g, eval_h }
00727 \endverbatim
00728 Note that the \c bool return flag for the evaluations methods 
00729 does not appear in the Ipopt documentation.
00730 Looking at the code, it seems to be a flag telling Ipopt to abort
00731 when the flag is false.
00732 
00733 */
00734 class cppad_ipopt_solution 
00735 {
00736 public:
00737      /// possible values for he solution status
00738      enum solution_status {
00739           not_defined,
00740           success,
00741           maxiter_exceeded,
00742           stop_at_tiny_step,
00743           stop_at_acceptable_point,
00744           local_infeasibility,
00745           user_requested_stop,
00746           feasible_point_found,
00747           diverging_iterates,
00748           restoration_failure,
00749           error_in_step_computation,
00750           invalid_number_detected,
00751           too_few_degrees_of_freedom,
00752           internal_error,
00753           unknown
00754      }  status;
00755      /// the approximation solution
00756      NumberVector      x;
00757      /// Lagrange multipliers corresponding to lower bounds on x
00758      NumberVector      z_l;
00759      /// Lagrange multipliers corresponding to upper bounds on x
00760      NumberVector      z_u;
00761      /// value of g(x)
00762      NumberVector      g;
00763      /// Lagrange multipliers correspondiing constraints on g(x)
00764      NumberVector      lambda;
00765      /// value of f(x)
00766      Ipopt::Number     obj_value;
00767      /// constructor initializes solution status as not yet defined
00768      cppad_ipopt_solution(void)
00769      {    status = not_defined; }
00770 };
00771 
00772 /*! 
00773 Class connects Ipopt to CppAD for derivative and sparsity pattern calculations.
00774 */
00775 class cppad_ipopt_nlp : public Ipopt::TNLP
00776 {
00777 private:
00778      /// A Scalar value used by Ipopt
00779      typedef Ipopt::Number                         Number;
00780      /// An index value used by Ipopt 
00781      typedef Ipopt::Index                          Index;
00782      /// Indexing style used in Ipopt sparsity structure 
00783      typedef Ipopt::TNLP::IndexStyleEnum           IndexStyleEnum;
00784      /// A simple vector of boolean values
00785      typedef CppAD::vectorBool                     BoolVector;
00786      /// A simple vector of AD function objects
00787      typedef CppAD::vector< CppAD::ADFun<Number> > ADFunVector;
00788      /// A simple vector of simple vectors of boolean values
00789      typedef CppAD::vector<BoolVector>             BoolVectorVector;
00790      /// A mapping that is dense in i, sparse in j, and maps (i, j)
00791      /// to the corresponding sparsity index in Ipopt.
00792      typedef CppAD::vector< std::map<size_t,size_t> > IndexMap;
00793 
00794      // ------------------------------------------------------------------
00795      // Values directly passed in to constuctor
00796      // ------------------------------------------------------------------
00797      /// dimension of the domain space for f(x) and g(x)
00798      /// (passed to ctor)
00799      const size_t                    n_;
00800      /// dimension of the range space for g(x)
00801      /// (passed to ctor)
00802      const size_t                    m_;
00803      /// dimension of the range space for g(x)
00804      /// (passed to ctor)
00805      const NumberVector              x_i_;
00806      /// lower limit for x 
00807      /// (size n_), (passed to ctor)
00808      const NumberVector              x_l_;
00809      /// upper limit for x 
00810      /// (size n_) (passed to ctor)
00811      const NumberVector              x_u_;
00812      /// lower limit for g(x)
00813      /// (size m_) (passed to ctor)
00814      const NumberVector              g_l_;
00815      /// upper limit for g(x)
00816      /// (size m_) (passed to ctor)
00817      const NumberVector              g_u_;
00818      /// pointer to base class version of derived class object used to get 
00819      /// information about the user's representation for f(x) and g(x)
00820      /// (passed to ctor)
00821      cppad_ipopt_fg_info* const      fg_info_;
00822      /// pointer to object where final results are stored
00823      /// (passed to ctor)
00824      cppad_ipopt_solution* const     solution_;
00825      /// plus infinity as a value of type Number
00826      const Number                    infinity_;
00827 
00828      // ------------------------------------------------------------------
00829      // Effectively const values determined during constructor using calls 
00830      // to fg_info:
00831      // ------------------------------------------------------------------
00832      /// The value of \f$ K \f$ in the representation.
00833      /// (effectively const)
00834      size_t                 K_;
00835      /// Does operation sequence for \f$ r_k (u) \f$ depend on \f$ u \f$.
00836      /// (size K_) (effectively const)
00837      BoolVector             retape_;
00838      /// <tt>q_[k]</tt> is the domain space dimension for \f$ r_k (u) \f$
00839      /// (size K_) (effectively const)
00840      SizeVector             q_;
00841      /// <tt>p_[k]</tt> is the range space dimension for \f$ r_k (u) \f$
00842      /// (size K_) (effectively const)
00843      SizeVector             p_;
00844      /// <tt>L_[k]</tt> is number of times \f$ r_k (u) \f$ appears in
00845      /// the representation summation
00846      /// (size K_) (effectively const)
00847      SizeVector             L_; 
00848      // -------------------------------------------------------------------
00849      // Other effectively const values determined by the constructor:
00850      // -------------------------------------------------------------------
00851      /*!
00852      CppAD sparsity patterns for \f$ \{ r_k^{(1)} (u) \} \f$ (set by ctor).
00853 
00854      For <tt>k = 0 , ... , K_-1, pattern_jac_r_[k]</tt> 
00855      is a CppAD sparsity pattern for the Jacobian of \f$ r_k (u) \f$ 
00856      and as such it has size <tt>p_[k]*q_[k]</tt>.
00857      (effectively const)
00858      */
00859      BoolVectorVector                 pattern_jac_r_;
00860 
00861      /*!
00862      CppAD sparsity patterns for \f$ \{ r_k^{(2)} (u) \} \f$ (set by ctor).
00863 
00864      For <tt>k = 0 , ... , K_-1, pattern_jac_r_[k]</tt> 
00865      is a CppAD sparsity pattern for the Hessian of 
00866      \f[
00867           R(u) = \sum_{i=0}^{p[k]-1}  r_k (u)_i
00868      \f]
00869      and as such it has size <tt>q_[k]*q_[k]</tt>.
00870      (effectively const)
00871      */
00872      BoolVectorVector                 pattern_hes_r_;
00873 
00874      /// number non-zero is Ipopt sparsity structor for Jacobian of g(x)
00875      /// (effectively const)
00876      size_t                           nnz_jac_g_;
00877      /// row indices in Ipopt sparsity structor for Jacobian of g(x)
00878      /// (effectively const)
00879      SizeVector                       iRow_jac_g_;
00880      /// column indices in Ipopt sparsity structor for Jacobian of g(x)
00881      /// (effectively const)
00882      SizeVector                       jCol_jac_g_;
00883 
00884      /// number non-zero is Ipopt sparsity structor for Hessian of Lagragian
00885      /// (effectively const)
00886      size_t                           nnz_h_lag_;
00887      /// row indices in Ipopt sparsity structor for Hessian of Lagragian
00888      /// (effectively const)
00889      SizeVector                       iRow_h_lag_;
00890      /// column indices in Ipopt sparsity structor for Hessian of Lagragian
00891      /// (effectively const)
00892      SizeVector                       jCol_h_lag_;
00893 
00894      /*!
00895      Mapping from (i, j) in Jacobian of g(x) to Ipopt sparsity structure
00896 
00897      For <tt>i = 0 , ... , m_-1, index_jac_g_[i]</tt>
00898      is a standard map from column index values \c j to the corresponding
00899      index in the Ipopt sparsity structure for the Jacobian of g(x). 
00900      */ 
00901      IndexMap                         index_jac_g_;
00902 
00903      /*!
00904      Mapping from (i, j) in Hessian of fg(x) to Ipopt sparsity structure
00905 
00906      For <tt>i = 0 , ... , n_-1, index_hes_fg_[i]</tt>
00907      is a standard map from column index values \c j to the corresponding
00908      index in the Ipopt sparsity structure for the Hessian of the Lagragian.
00909      */ 
00910      IndexMap                         index_hes_fg_;
00911      // -----------------------------------------------------------------
00912      // Values that are changed by routine other than the constructor:
00913      // -----------------------------------------------------------------
00914 
00915      /// For <tt>k = 0 , ... , K_-1, r_fun_[k]</tt> 
00916      /// is a the CppAD function object corresponding to \f$ r_k (u) \f$.
00917      ADFunVector                      r_fun_;
00918      /*!
00919      Is r_fun[k] OK for current x.
00920  
00921      For <tt>k = 0 , ... , K_-1, tape_ok_[k]</tt> 
00922      is true if current operations sequence in <tt>r_fun_[k]</tt> 
00923      OK for this value of \f$ x \f$. 
00924      Note that \f$ u = [ J_{k,\ell} \otimes n ] (x) \f$ may depend on the
00925      value of \f$ \ell \f$.
00926      */
00927      BoolVector             tape_ok_;
00928 
00929      /// work space of size equal maximum of <tt>q[k]</tt> w.r.t \c k.
00930      SizeVector             J_;
00931      /// work space of size equal maximum of <tt>p[k]</tt> w.r.t \c k.
00932      SizeVector             I_; 
00933      // ------------------------------------------------------------
00934      // Private Methods
00935      // ------------------------------------------------------------
00936      /// block the default constructor from use
00937      cppad_ipopt_nlp(const cppad_ipopt_nlp&);
00938      /// blocks the assignment operator from use
00939      cppad_ipopt_nlp& operator=(const cppad_ipopt_nlp&);
00940 public:
00941      // ----------------------------------------------------------------
00942      // See cppad_ipopt_nlp.cpp for doxygen documentation of these methods
00943      // ----------------------------------------------------------------
00944 
00945      /// only constructor for cppad_ipopot_nlp
00946      cppad_ipopt_nlp(
00947           size_t n                         , 
00948           size_t m                         ,
00949           const NumberVector    &x_i       ,
00950           const NumberVector    &x_l       ,
00951           const NumberVector    &x_u       ,
00952           const NumberVector    &g_l       ,
00953           const NumberVector    &g_u       ,
00954           cppad_ipopt_fg_info*   fg_info   ,
00955           cppad_ipopt_solution*  solution
00956      );
00957 
00958      // use virtual so that derived class destructor gets called.
00959      virtual ~cppad_ipopt_nlp();
00960 
00961      // return info about the nlp
00962      virtual bool get_nlp_info(
00963           Index&          n           , 
00964           Index&          m           , 
00965           Index&          nnz_jac_g   ,
00966           Index&          nnz_h_lag   , 
00967           IndexStyleEnum& index_style
00968      );
00969 
00970      // return bounds for my problem 
00971      virtual bool get_bounds_info(
00972           Index           n   , 
00973           Number*         x_l , 
00974           Number*         x_u ,
00975           Index           m   , 
00976           Number*         g_l , 
00977           Number*         g_u   
00978      );
00979 
00980      // return the starting point for the algorithm 
00981      virtual bool get_starting_point(
00982           Index          n            , 
00983           bool           init_x       , 
00984           Number*        x            ,
00985           bool           init_z       , 
00986           Number*        z_L          , 
00987           Number*        z_U          ,
00988           Index          m            ,
00989           bool           init_lambda  ,
00990           Number*        lambda
00991      );
00992 
00993      // return the objective value 
00994      virtual bool eval_f(
00995           Index          n           , 
00996           const Number*  x           , 
00997           bool           new_x       , 
00998           Number&        obj_value
00999      );
01000 
01001      // Method to return the gradient of the objective 
01002      virtual bool eval_grad_f(
01003           Index          n           , 
01004           const Number*  x           , 
01005           bool           new_x       , 
01006           Number*        grad_f
01007      );
01008 
01009      // return the constraint residuals
01010      virtual bool eval_g(
01011           Index          n           , 
01012           const Number*  x           , 
01013           bool           new_x       , 
01014           Index          m           , 
01015           Number*        g
01016      );
01017 
01018      // Method to return:
01019      // 1) The structure of the jacobian (if "values" is NULL)
01020      // 2) The values of the jacobian (if "values" is not NULL)
01021      virtual bool eval_jac_g(
01022           Index          n           , 
01023           const Number*  x           , 
01024           bool           new_x       ,
01025           Index          m           , 
01026           Index          nele_jac    , 
01027           Index*         iRow        , 
01028           Index*         jCol        ,
01029           Number*        values
01030      );
01031 
01032      // Method to return:
01033      //  1) structure of hessian of the lagrangian (if "values" is NULL)
01034      //  2) values of hessian of the lagrangian (if "values" is not NULL)
01035      virtual bool eval_h(
01036           Index          n           , 
01037           const Number*  x           , 
01038           bool           new_x       ,
01039           Number         obj_factor  , 
01040           Index          m           , 
01041           const Number*  lambda      ,
01042           bool           new_lambda  , 
01043           Index          nele_hess   , 
01044           Index*         iRow        ,
01045           Index*         jCol        , 
01046           Number*        values
01047      );
01048 
01049      // called when the algorithm is completed so the TNLP can 
01050      // store/write the solution 
01051      virtual void finalize_solution(
01052           Ipopt::SolverReturn       status      ,
01053           Index                      n          , 
01054           const Number*              x          , 
01055           const Number*              z_L        , 
01056           const Number*              z_U        ,
01057           Index                      m          , 
01058           const Number*              g          , 
01059           const Number*              lambda     ,
01060           Number                     obj_value  ,
01061           const Ipopt::IpoptData*           ip_data    ,
01062           Ipopt::IpoptCalculatedQuantities* ip_cq
01063      );
01064 
01065      virtual bool intermediate_callback(
01066           Ipopt::AlgorithmMode              mode,
01067           Index                             iter, 
01068           Number                            obj_value,
01069           Number                            inf_pr, 
01070           Number                            inf_du,
01071           Number                            mu, 
01072           Number                            d_norm,
01073           Number                            regularization_size,
01074           Number                            alpha_du, 
01075           Number                            alpha_pr,
01076           Index                             ls_trials,
01077           const Ipopt::IpoptData*           ip_data,
01078           Ipopt::IpoptCalculatedQuantities* ip_cq
01079      );
01080 
01081 };
01082 
01083 
01084 // ---------------------------------------------------------------------------
01085 } // end namespace cppad_ipopt
01086 // ---------------------------------------------------------------------------
01087 
01088 # endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines