00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #ifndef __MITTELMANNBNDRYCNTRLNEUM_HPP__
00011 #define __MITTELMANNBNDRYCNTRLNEUM_HPP__
00012
00013 #include "IpTNLP.hpp"
00014 #include "RegisteredTNLP.hpp"
00015
00016 #ifdef HAVE_CMATH
00017 # include <cmath>
00018 #else
00019 # ifdef HAVE_MATH_H
00020 # include <math.h>
00021 # else
00022 # error "don't have header file for math"
00023 # endif
00024 #endif
00025
00026 using namespace Ipopt;
00027
00033 class MittelmannBndryCntrlNeumBase : public RegisteredTNLP
00034 {
00035 public:
00038 MittelmannBndryCntrlNeumBase();
00039
00041 virtual ~MittelmannBndryCntrlNeumBase();
00042
00046 virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00047 Index& nnz_h_lag, IndexStyleEnum& index_style);
00048
00050 virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00051 Index m, Number* g_l, Number* g_u);
00052
00054 virtual bool get_starting_point(Index n, bool init_x, Number* x,
00055 bool init_z, Number* z_L, Number* z_U,
00056 Index m, bool init_lambda,
00057 Number* lambda);
00058
00060 virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00061
00063 virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00064
00066 virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00067
00072 virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00073 Index m, Index nele_jac, Index* iRow, Index *jCol,
00074 Number* values);
00075
00080 virtual bool eval_h(Index n, const Number* x, bool new_x,
00081 Number obj_factor, Index m, const Number* lambda,
00082 bool new_lambda, Index nele_hess, Index* iRow,
00083 Index* jCol, Number* values);
00084
00086
00088 virtual bool get_scaling_parameters(Number& obj_scaling,
00089 bool& use_x_scaling, Index n,
00090 Number* x_scaling,
00091 bool& use_g_scaling, Index m,
00092 Number* g_scaling);
00093
00098 virtual void finalize_solution(SolverReturn status,
00099 Index n, const Number* x, const Number* z_L, const Number* z_U,
00100 Index m, const Number* g, const Number* lambda,
00101 Number obj_value,
00102 const IpoptData* ip_data,
00103 IpoptCalculatedQuantities* ip_cq);
00105
00106 protected:
00110 void SetBaseParameters(Index N, Number alpha, Number lb_y,
00111 Number ub_y, Number lb_u, Number ub_u,
00112 Number u_init);
00113
00117 virtual Number y_d_cont(Number x1, Number x2) const =0;
00119 virtual Number d_cont(Number x1, Number x2, Number y) const =0;
00121 virtual Number d_cont_dy(Number x1, Number x2, Number y) const =0;
00123 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const =0;
00126 virtual bool d_cont_dydy_alwayszero() const =0;
00128 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const =0;
00130 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00132 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00134 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00137 virtual bool b_cont_dydy_alwayszero() const =0;
00139
00140 private:
00152 MittelmannBndryCntrlNeumBase(const MittelmannBndryCntrlNeumBase&);
00153 MittelmannBndryCntrlNeumBase& operator=(const MittelmannBndryCntrlNeumBase&);
00155
00159 Index N_;
00161 Number h_;
00163 Number hh_;
00165 Number lb_y_;
00167 Number ub_y_;
00169 Number lb_u_;
00171 Number ub_u_;
00173 Number u_init_;
00176 Number alpha_;
00178 Number* y_d_;
00180
00185 inline Index y_index(Index i, Index j) const
00186 {
00187 return j + (N_+2)*i;
00188 }
00191 inline Index u0j_index(Index j) const
00192 {
00193 return (N_+2)*(N_+2) + j-1;
00194 }
00197 inline Index u1j_index(Index j) const
00198 {
00199 return (N_+2)*(N_+2) + N_ + j-1;
00200 }
00203 inline Index ui0_index(Index j) const
00204 {
00205 return (N_+2)*(N_+2) + 2*N_ + j-1;
00206 }
00209 inline Index ui1_index(Index j) const
00210 {
00211 return (N_+2)*(N_+2) + 3*N_ + j-1;
00212 }
00214 inline Number x1_grid(Index i) const
00215 {
00216 return h_*(Number)i;
00217 }
00219 inline Number x2_grid(Index j) const
00220 {
00221 return h_*(Number)j;
00222 }
00224 };
00225
00227 class MittelmannBndryCntrlNeum1 : public MittelmannBndryCntrlNeumBase
00228 {
00229 public:
00230 MittelmannBndryCntrlNeum1()
00231 {}
00232
00233 virtual ~MittelmannBndryCntrlNeum1()
00234 {}
00235
00236 virtual bool InitializeProblem(Index N)
00237 {
00238 if (N<1) {
00239 printf("N has to be at least 1.");
00240 return false;
00241 }
00242 Number alpha = 0.01;
00243 Number lb_y = -1e20;
00244 Number ub_y = 2.071;
00245 Number lb_u = 3.7;
00246 Number ub_u = 4.5;
00247 Number u_init = (ub_u+lb_u)/2.;
00248
00249 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00250 return true;
00251 }
00252 protected:
00254 virtual Number y_d_cont(Number x1, Number x2) const
00255 {
00256 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00257 }
00259 virtual Number d_cont(Number x1, Number x2, Number y) const
00260 {
00261 return 0.;
00262 }
00264 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00265 {
00266 return 0.;
00267 }
00269 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00270 {
00271 return 0.;
00272 }
00275 virtual bool d_cont_dydy_alwayszero() const
00276 {
00277 return true;
00278 }
00280 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00281 {
00282 return u - y*y;
00283 }
00285 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00286 {
00287 return - 2.*y;
00288 }
00290 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00291 {
00292 return 1.;
00293 }
00295 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00296 {
00297 return -2.;
00298 }
00301 virtual bool b_cont_dydy_alwayszero() const
00302 {
00303 return false;
00304 }
00305 private:
00308 MittelmannBndryCntrlNeum1(const MittelmannBndryCntrlNeum1&);
00309 MittelmannBndryCntrlNeum1& operator=(const MittelmannBndryCntrlNeum1&);
00311 };
00312
00314 class MittelmannBndryCntrlNeum2 : public MittelmannBndryCntrlNeumBase
00315 {
00316 public:
00317 MittelmannBndryCntrlNeum2()
00318 {}
00319
00320 virtual ~MittelmannBndryCntrlNeum2()
00321 {}
00322
00323 virtual bool InitializeProblem(Index N)
00324 {
00325 if (N<1) {
00326 printf("N has to be at least 1.");
00327 return false;
00328 }
00329 Number alpha = 0.;
00330 Number lb_y = -1e20;
00331 Number ub_y = 2.835;
00332 Number lb_u = 6.;
00333 Number ub_u = 9.;
00334 Number u_init = (ub_u+lb_u)/2.;
00335
00336 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00337 return true;
00338 }
00339 protected:
00341 virtual Number y_d_cont(Number x1, Number x2) const
00342 {
00343 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00344 }
00346 virtual Number d_cont(Number x1, Number x2, Number y) const
00347 {
00348 return 0.;
00349 }
00351 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00352 {
00353 return 0.;
00354 }
00356 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00357 {
00358 return 0.;
00359 }
00362 virtual bool d_cont_dydy_alwayszero() const
00363 {
00364 return true;
00365 }
00367 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00368 {
00369 return u - y*y;
00370 }
00372 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00373 {
00374 return - 2.*y;
00375 }
00377 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00378 {
00379 return 1.;
00380 }
00382 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00383 {
00384 return -2.;
00385 }
00388 virtual bool b_cont_dydy_alwayszero() const
00389 {
00390 return false;
00391 }
00392 private:
00395 MittelmannBndryCntrlNeum2(const MittelmannBndryCntrlNeum2&);
00396 MittelmannBndryCntrlNeum2& operator=(const MittelmannBndryCntrlNeum2&);
00398 };
00399
00401 class MittelmannBndryCntrlNeum3 : public MittelmannBndryCntrlNeumBase
00402 {
00403 public:
00404 MittelmannBndryCntrlNeum3()
00405 {}
00406
00407 virtual ~MittelmannBndryCntrlNeum3()
00408 {}
00409
00410 virtual bool InitializeProblem(Index N)
00411 {
00412 if (N<1) {
00413 printf("N has to be at least 1.");
00414 return false;
00415 }
00416 Number alpha = 0.01;
00417 Number lb_y = -1e20;
00418 Number ub_y = 2.7;
00419 Number lb_u = 1.8;
00420 Number ub_u = 2.5;
00421 Number u_init = (ub_u+lb_u)/2.;
00422
00423 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00424 return true;
00425 }
00426 protected:
00428 virtual Number y_d_cont(Number x1, Number x2) const
00429 {
00430 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00431 }
00433 virtual Number d_cont(Number x1, Number x2, Number y) const
00434 {
00435 return y*y*y-y;
00436 }
00438 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00439 {
00440 return 3.*y*y-1.;
00441 }
00443 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00444 {
00445 return 6.*y;
00446 }
00449 virtual bool d_cont_dydy_alwayszero() const
00450 {
00451 return false;
00452 }
00454 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00455 {
00456 return u;
00457 }
00459 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00460 {
00461 return 0.;
00462 }
00464 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00465 {
00466 return 1.;
00467 }
00469 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00470 {
00471 return 0.;
00472 }
00475 virtual bool b_cont_dydy_alwayszero() const
00476 {
00477 return true;
00478 }
00479 private:
00482 MittelmannBndryCntrlNeum3(const MittelmannBndryCntrlNeum3&);
00483 MittelmannBndryCntrlNeum3& operator=(const MittelmannBndryCntrlNeum3&);
00485 };
00486
00488 class MittelmannBndryCntrlNeum4 : public MittelmannBndryCntrlNeumBase
00489 {
00490 public:
00491 MittelmannBndryCntrlNeum4()
00492 {}
00493
00494 virtual ~MittelmannBndryCntrlNeum4()
00495 {}
00496
00497 virtual bool InitializeProblem(Index N)
00498 {
00499 if (N<1) {
00500 printf("N has to be at least 1.");
00501 return false;
00502 }
00503 Number alpha = 0.;
00504 Number lb_y = -1e20;
00505 Number ub_y = 2.7;
00506 Number lb_u = 1.8;
00507 Number ub_u = 2.5;
00508 Number u_init = (ub_u+lb_u)/2.;
00509
00510 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00511 return true;
00512 }
00513 protected:
00515 virtual Number y_d_cont(Number x1, Number x2) const
00516 {
00517 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00518 }
00520 virtual Number d_cont(Number x1, Number x2, Number y) const
00521 {
00522 return y*y*y-y;
00523 }
00525 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00526 {
00527 return 3.*y*y-1.;
00528 }
00530 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00531 {
00532 return 6.*y;
00533 }
00536 virtual bool d_cont_dydy_alwayszero() const
00537 {
00538 return false;
00539 }
00541 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00542 {
00543 return u;
00544 }
00546 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00547 {
00548 return 0.;
00549 }
00551 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00552 {
00553 return 1.;
00554 }
00556 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00557 {
00558 return 0.;
00559 }
00562 virtual bool b_cont_dydy_alwayszero() const
00563 {
00564 return true;
00565 }
00566 private:
00569 MittelmannBndryCntrlNeum4(const MittelmannBndryCntrlNeum4&);
00570 MittelmannBndryCntrlNeum4& operator=(const MittelmannBndryCntrlNeum4&);
00572 };
00573
00574 #endif