Ipopt  trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
MittelmannDistCntrlNeumA.hpp
Go to the documentation of this file.
00001 // Copyright (C) 2005, 2006 International Business Machines and others.
00002 // All Rights Reserved.
00003 // This code is published under the Eclipse Public License.
00004 //
00005 // $Id$
00006 //
00007 // Authors:  Andreas Waechter           IBM    2005-10-18
00008 //                  based on MyNLP.hpp
00009 
00010 #ifndef __MITTELMANNDISTRCNTRLNEUMA_HPP__
00011 #define __MITTELMANNDISTRCNTRLNEUMA_HPP__
00012 
00013 #include "IpTNLP.hpp"
00014 #include "RegisteredTNLP.hpp"
00015 
00016 #ifdef HAVE_CONFIG_H
00017 #include "config.h"
00018 #else
00019 #include "configall_system.h"
00020 #endif
00021 
00022 #ifdef HAVE_CMATH
00023 # include <cmath>
00024 #else
00025 # ifdef HAVE_MATH_H
00026 #  include <math.h>
00027 # else
00028 #  error "don't have header file for math"
00029 # endif
00030 #endif
00031 
00032 #ifdef HAVE_CSTDIO
00033 # include <cstdio>
00034 #else
00035 # ifdef HAVE_STDIO_H
00036 #  include <stdio.h>
00037 # else
00038 #  error "don't have header file for stdio"
00039 # endif
00040 #endif
00041 
00042 using namespace Ipopt;
00043 
00050 class MittelmannDistCntrlNeumABase : public RegisteredTNLP
00051 {
00052 public:
00055   MittelmannDistCntrlNeumABase();
00056 
00058   virtual ~MittelmannDistCntrlNeumABase();
00059 
00063   virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00064                             Index& nnz_h_lag, IndexStyleEnum& index_style);
00065 
00067   virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00068                                Index m, Number* g_l, Number* g_u);
00069 
00071   virtual bool get_starting_point(Index n, bool init_x, Number* x,
00072                                   bool init_z, Number* z_L, Number* z_U,
00073                                   Index m, bool init_lambda,
00074                                   Number* lambda);
00075 
00077   virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00078 
00080   virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00081 
00083   virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00084 
00089   virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00090                           Index m, Index nele_jac, Index* iRow, Index *jCol,
00091                           Number* values);
00092 
00097   virtual bool eval_h(Index n, const Number* x, bool new_x,
00098                       Number obj_factor, Index m, const Number* lambda,
00099                       bool new_lambda, Index nele_hess, Index* iRow,
00100                       Index* jCol, Number* values);
00101 
00103 
00105   virtual bool get_scaling_parameters(Number& obj_scaling,
00106                                       bool& use_x_scaling, Index n,
00107                                       Number* x_scaling,
00108                                       bool& use_g_scaling, Index m,
00109                                       Number* g_scaling);
00110 
00115   virtual void finalize_solution(SolverReturn status,
00116                                  Index n, const Number* x, const Number* z_L, const Number* z_U,
00117                                  Index m, const Number* g, const Number* lambda,
00118                                  Number obj_value,
00119                                  const IpoptData* ip_data,
00120                                  IpoptCalculatedQuantities* ip_cq);
00122 
00123 protected:
00127   void SetBaseParameters(Index N, Number lb_y,
00128                          Number ub_y, Number lb_u, Number ub_u,
00129                          Number b_0j, Number b_1j, Number b_i0, Number b_i1,
00130                          Number u_init);
00131 
00135   virtual Number y_d_cont(Number x1, Number x2) const =0;
00137   virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const =0;
00139   virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00141   virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00143   virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00146   virtual bool fint_cont_dydy_alwayszero() const =0;
00148   virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00151   virtual bool fint_cont_dudu_alwayszero() const =0;
00153   virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00156   virtual bool fint_cont_dydu_alwayszero() const =0;
00158   virtual Number d_cont(Number x1, Number x2, Number y, Number u) const =0;
00160   virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00162   virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00164   virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00167   virtual bool d_cont_dydy_alwayszero() const =0;
00169   virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00172   virtual bool d_cont_dudu_alwayszero() const =0;
00174   virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00177   virtual bool d_cont_dydu_alwayszero() const =0;
00179 
00180 private:
00192   MittelmannDistCntrlNeumABase(const MittelmannDistCntrlNeumABase&);
00193   MittelmannDistCntrlNeumABase& operator=(const MittelmannDistCntrlNeumABase&);
00195 
00199   Index N_;
00201   Number h_;
00203   Number hh_;
00205   Number lb_y_;
00207   Number ub_y_;
00209   Number lb_u_;
00211   Number ub_u_;
00214   Number b_0j_;
00217   Number b_1j_;
00220   Number b_i0_;
00223   Number b_i1_;
00225   Number u_init_;
00227   Number* y_d_;
00229 
00234   inline Index y_index(Index i, Index j) const
00235   {
00236     return j + (N_+2)*i;
00237   }
00240   inline Index u_index(Index i, Index j) const
00241   {
00242     return (N_+2)*(N_+2) + (j-1) + (N_)*(i-1);
00243   }
00246   inline Index pde_index(Index i, Index j) const
00247   {
00248     return (j-1) + N_*(i-1);
00249   }
00251   inline Number x1_grid(Index i) const
00252   {
00253     return h_*(Number)i;
00254   }
00256   inline Number x2_grid(Index i) const
00257   {
00258     return h_*(Number)i;
00259   }
00261 };
00262 
00264 class MittelmannDistCntrlNeumA1 : public MittelmannDistCntrlNeumABase
00265 {
00266 public:
00267   MittelmannDistCntrlNeumA1()
00268       :
00269       pi_(4.*atan(1.)),
00270       alpha_(0.001)
00271   {}
00272 
00273   virtual ~MittelmannDistCntrlNeumA1()
00274   {}
00275 
00276   virtual bool InitializeProblem(Index N)
00277   {
00278     if (N<1) {
00279       printf("N has to be at least 1.");
00280       return false;
00281     }
00282     Number lb_y = -1e20;
00283     Number ub_y = 0.371;
00284     Number lb_u = -8.;
00285     Number ub_u = 9.;
00286     Number b_0j = 1.;
00287     Number b_1j = 1.;
00288     Number b_i0 = 1.;
00289     Number b_i1 = 1.;
00290     Number u_init = (ub_u+lb_u)/2.;
00291 
00292     SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00293     return true;
00294   }
00295 protected:
00297   virtual Number y_d_cont(Number x1, Number x2)  const
00298   {
00299     return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00300   }
00302   virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00303   {
00304     Number diff_y = y-y_d_cont(x1,x2);
00305     return 0.5*(diff_y*diff_y + alpha_*u*u);
00306   }
00308   virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00309   {
00310     return  y-y_d_cont(x1,x2);
00311   }
00312 
00314   virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00315   {
00316     return alpha_*u;
00317   }
00319   virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00320   {
00321     return 1.;
00322   }
00325   virtual bool fint_cont_dydy_alwayszero() const
00326   {
00327     return false;
00328   }
00330   virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00331   {
00332     return alpha_;
00333   }
00336   virtual bool fint_cont_dudu_alwayszero() const
00337   {
00338     return false;
00339   }
00341   virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00342   {
00343     return 0.;
00344   }
00347   virtual bool fint_cont_dydu_alwayszero() const
00348   {
00349     return true;
00350   }
00352   virtual Number d_cont(Number x1, Number x2, Number y, Number u)  const
00353   {
00354     return -exp(y) - u;
00355   }
00357   virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u)  const
00358   {
00359     return -exp(y);
00360   }
00362   virtual Number d_cont_du(Number x1, Number x2, Number y, Number u)  const
00363   {
00364     return -1.;
00365   }
00367   virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u)  const
00368   {
00369     return -exp(y);
00370   }
00373   virtual bool d_cont_dydy_alwayszero() const
00374   {
00375     return false;
00376   }
00378   virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00379   {
00380     return 0.;
00381   }
00384   virtual bool d_cont_dudu_alwayszero() const
00385   {
00386     return true;
00387   }
00389   virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00390   {
00391     return 0.;
00392   }
00395   virtual bool d_cont_dydu_alwayszero() const
00396   {
00397     return true;
00398   }
00399 private:
00402   MittelmannDistCntrlNeumA1(const MittelmannDistCntrlNeumA1&);
00403   MittelmannDistCntrlNeumA1& operator=(const MittelmannDistCntrlNeumA1&);
00405 
00406   const Number pi_;
00408   const Number alpha_;
00409 };
00410 
00412 class MittelmannDistCntrlNeumA2 : public MittelmannDistCntrlNeumABase
00413 {
00414 public:
00415   MittelmannDistCntrlNeumA2()
00416       :
00417       pi_(4.*atan(1.))
00418   {}
00419 
00420   virtual ~MittelmannDistCntrlNeumA2()
00421   {}
00422 
00423   virtual bool InitializeProblem(Index N)
00424   {
00425     if (N<1) {
00426       printf("N has to be at least 1.");
00427       return false;
00428     }
00429     Number lb_y = -1e20;
00430     Number ub_y = 0.371;
00431     Number lb_u = -8.;
00432     Number ub_u = 9.;
00433     Number b_0j = 1.;
00434     Number b_1j = 1.;
00435     Number b_i0 = 1.;
00436     Number b_i1 = 1.;
00437     Number u_init = (ub_u+lb_u)/2.;
00438 
00439     SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00440     return true;
00441   }
00442 protected:
00444   virtual Number y_d_cont(Number x1, Number x2)  const
00445   {
00446     return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00447   }
00449   virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00450   {
00451     Number diff_y = y-y_d_cont(x1,x2);
00452     return 0.5*diff_y*diff_y;
00453   }
00455   virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00456   {
00457     return  y-y_d_cont(x1,x2);
00458   }
00459 
00461   virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00462   {
00463     return 0.;
00464   }
00466   virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00467   {
00468     return 1.;
00469   }
00472   virtual bool fint_cont_dydy_alwayszero() const
00473   {
00474     return false;
00475   }
00477   virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00478   {
00479     return 0.;
00480   }
00483   virtual bool fint_cont_dudu_alwayszero() const
00484   {
00485     return true;
00486   }
00488   virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00489   {
00490     return 0.;
00491   }
00494   virtual bool fint_cont_dydu_alwayszero() const
00495   {
00496     return true;
00497   }
00499   virtual Number d_cont(Number x1, Number x2, Number y, Number u)  const
00500   {
00501     return -exp(y) - u;
00502   }
00504   virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u)  const
00505   {
00506     return -exp(y);
00507   }
00509   virtual Number d_cont_du(Number x1, Number x2, Number y, Number u)  const
00510   {
00511     return -1.;
00512   }
00514   virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u)  const
00515   {
00516     return -exp(y);
00517   }
00520   virtual bool d_cont_dydy_alwayszero() const
00521   {
00522     return false;
00523   }
00525   virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00526   {
00527     return 0.;
00528   }
00531   virtual bool d_cont_dudu_alwayszero() const
00532   {
00533     return true;
00534   }
00536   virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00537   {
00538     return 0.;
00539   }
00542   virtual bool d_cont_dydu_alwayszero() const
00543   {
00544     return true;
00545   }
00546 private:
00549   MittelmannDistCntrlNeumA2(const MittelmannDistCntrlNeumA2&);
00550   MittelmannDistCntrlNeumA2& operator=(const MittelmannDistCntrlNeumA2&);
00552 
00553   const Number pi_;
00554 };
00555 
00557 class MittelmannDistCntrlNeumA3 : public MittelmannDistCntrlNeumABase
00558 {
00559 public:
00560   MittelmannDistCntrlNeumA3()
00561       :
00562       pi_(4.*atan(1.)),
00563       M_(1.),
00564       K_(0.8),
00565       b_(1.)
00566   {}
00567 
00568   virtual ~MittelmannDistCntrlNeumA3()
00569   {}
00570 
00571   virtual bool InitializeProblem(Index N)
00572   {
00573     if (N<1) {
00574       printf("N has to be at least 1.");
00575       return false;
00576     }
00577     Number lb_y = 3.;//-1e20;
00578     Number ub_y = 6.09;
00579     Number lb_u = 1.4;
00580     Number ub_u = 1.6;
00581     Number b_0j = 1.;
00582     Number b_1j = 0.;
00583     Number b_i0 = 1.;
00584     Number b_i1 = 0.;
00585     Number u_init = (ub_u+lb_u)/2.;
00586 
00587     SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00588     return true;
00589   }
00590 protected:
00592   virtual Number y_d_cont(Number x1, Number x2)  const
00593   {
00594     return 6.;
00595   }
00597   virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00598   {
00599     return u*(M_*u - K_*y);
00600   }
00602   virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00603   {
00604     return -K_*u;
00605   }
00606 
00608   virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00609   {
00610     return 2.*M_*u - K_*y;
00611   }
00613   virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00614   {
00615     return 0.;
00616   }
00619   virtual bool fint_cont_dydy_alwayszero() const
00620   {
00621     return true;
00622   }
00624   virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00625   {
00626     return 2.*M_;
00627   }
00630   virtual bool fint_cont_dudu_alwayszero() const
00631   {
00632     return false;
00633   }
00635   virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00636   {
00637     return -K_;
00638   }
00641   virtual bool fint_cont_dydu_alwayszero() const
00642   {
00643     return false;
00644   }
00646   virtual Number d_cont(Number x1, Number x2, Number y, Number u)  const
00647   {
00648     return y*(u + b_*y - a(x1,x2));
00649   }
00651   virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u)  const
00652   {
00653     return (u + 2.*b_*y -a(x1,x2));
00654   }
00656   virtual Number d_cont_du(Number x1, Number x2, Number y, Number u)  const
00657   {
00658     return y;
00659   }
00661   virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u)  const
00662   {
00663     return 2.*b_;
00664   }
00667   virtual bool d_cont_dydy_alwayszero() const
00668   {
00669     return false;
00670   }
00672   virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00673   {
00674     return 0.;
00675   }
00678   virtual bool d_cont_dudu_alwayszero() const
00679   {
00680     return true;
00681   }
00683   virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00684   {
00685     return 1.;
00686   }
00689   virtual bool d_cont_dydu_alwayszero() const
00690   {
00691     return false;
00692   }
00693 private:
00696   MittelmannDistCntrlNeumA3(const MittelmannDistCntrlNeumA3&);
00697   MittelmannDistCntrlNeumA3& operator=(const MittelmannDistCntrlNeumA3&);
00699 
00700   const Number pi_;
00701   /*@name constrants appearing in problem formulation */
00703   const Number M_;
00704   const Number K_;
00705   const Number b_;
00707   //* Auxiliary function for state equation */
00708   inline Number a(Number x1, Number x2) const
00709   {
00710     return 7. + 4.*sin(2.*pi_*x1*x2);
00711   }
00712 };
00713 
00714 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines