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_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
00049 class MittelmannBndryCntrlNeumBase : public RegisteredTNLP
00050 {
00051 public:
00054 MittelmannBndryCntrlNeumBase();
00055
00057 virtual ~MittelmannBndryCntrlNeumBase();
00058
00062 virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00063 Index& nnz_h_lag, IndexStyleEnum& index_style);
00064
00066 virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00067 Index m, Number* g_l, Number* g_u);
00068
00070 virtual bool get_starting_point(Index n, bool init_x, Number* x,
00071 bool init_z, Number* z_L, Number* z_U,
00072 Index m, bool init_lambda,
00073 Number* lambda);
00074
00076 virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00077
00079 virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00080
00082 virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00083
00088 virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00089 Index m, Index nele_jac, Index* iRow, Index *jCol,
00090 Number* values);
00091
00096 virtual bool eval_h(Index n, const Number* x, bool new_x,
00097 Number obj_factor, Index m, const Number* lambda,
00098 bool new_lambda, Index nele_hess, Index* iRow,
00099 Index* jCol, Number* values);
00100
00102
00104 virtual bool get_scaling_parameters(Number& obj_scaling,
00105 bool& use_x_scaling, Index n,
00106 Number* x_scaling,
00107 bool& use_g_scaling, Index m,
00108 Number* g_scaling);
00109
00114 virtual void finalize_solution(SolverReturn status,
00115 Index n, const Number* x, const Number* z_L, const Number* z_U,
00116 Index m, const Number* g, const Number* lambda,
00117 Number obj_value,
00118 const IpoptData* ip_data,
00119 IpoptCalculatedQuantities* ip_cq);
00121
00122 protected:
00126 void SetBaseParameters(Index N, Number alpha, Number lb_y,
00127 Number ub_y, Number lb_u, Number ub_u,
00128 Number u_init);
00129
00133 virtual Number y_d_cont(Number x1, Number x2) const =0;
00135 virtual Number d_cont(Number x1, Number x2, Number y) const =0;
00137 virtual Number d_cont_dy(Number x1, Number x2, Number y) const =0;
00139 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const =0;
00142 virtual bool d_cont_dydy_alwayszero() const =0;
00144 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const =0;
00146 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00148 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00150 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00153 virtual bool b_cont_dydy_alwayszero() const =0;
00155
00156 private:
00168 MittelmannBndryCntrlNeumBase(const MittelmannBndryCntrlNeumBase&);
00169 MittelmannBndryCntrlNeumBase& operator=(const MittelmannBndryCntrlNeumBase&);
00171
00175 Index N_;
00177 Number h_;
00179 Number hh_;
00181 Number lb_y_;
00183 Number ub_y_;
00185 Number lb_u_;
00187 Number ub_u_;
00189 Number u_init_;
00192 Number alpha_;
00194 Number* y_d_;
00196
00201 inline Index y_index(Index i, Index j) const
00202 {
00203 return j + (N_+2)*i;
00204 }
00207 inline Index u0j_index(Index j) const
00208 {
00209 return (N_+2)*(N_+2) + j-1;
00210 }
00213 inline Index u1j_index(Index j) const
00214 {
00215 return (N_+2)*(N_+2) + N_ + j-1;
00216 }
00219 inline Index ui0_index(Index j) const
00220 {
00221 return (N_+2)*(N_+2) + 2*N_ + j-1;
00222 }
00225 inline Index ui1_index(Index j) const
00226 {
00227 return (N_+2)*(N_+2) + 3*N_ + j-1;
00228 }
00230 inline Number x1_grid(Index i) const
00231 {
00232 return h_*(Number)i;
00233 }
00235 inline Number x2_grid(Index j) const
00236 {
00237 return h_*(Number)j;
00238 }
00240 };
00241
00243 class MittelmannBndryCntrlNeum1 : public MittelmannBndryCntrlNeumBase
00244 {
00245 public:
00246 MittelmannBndryCntrlNeum1()
00247 {}
00248
00249 virtual ~MittelmannBndryCntrlNeum1()
00250 {}
00251
00252 virtual bool InitializeProblem(Index N)
00253 {
00254 if (N<1) {
00255 printf("N has to be at least 1.");
00256 return false;
00257 }
00258 Number alpha = 0.01;
00259 Number lb_y = -1e20;
00260 Number ub_y = 2.071;
00261 Number lb_u = 3.7;
00262 Number ub_u = 4.5;
00263 Number u_init = (ub_u+lb_u)/2.;
00264
00265 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00266 return true;
00267 }
00268 protected:
00270 virtual Number y_d_cont(Number x1, Number x2) const
00271 {
00272 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00273 }
00275 virtual Number d_cont(Number x1, Number x2, Number y) const
00276 {
00277 return 0.;
00278 }
00280 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00281 {
00282 return 0.;
00283 }
00285 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00286 {
00287 return 0.;
00288 }
00291 virtual bool d_cont_dydy_alwayszero() const
00292 {
00293 return true;
00294 }
00296 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00297 {
00298 return u - y*y;
00299 }
00301 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00302 {
00303 return - 2.*y;
00304 }
00306 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00307 {
00308 return 1.;
00309 }
00311 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00312 {
00313 return -2.;
00314 }
00317 virtual bool b_cont_dydy_alwayszero() const
00318 {
00319 return false;
00320 }
00321 private:
00324 MittelmannBndryCntrlNeum1(const MittelmannBndryCntrlNeum1&);
00325 MittelmannBndryCntrlNeum1& operator=(const MittelmannBndryCntrlNeum1&);
00327 };
00328
00330 class MittelmannBndryCntrlNeum2 : public MittelmannBndryCntrlNeumBase
00331 {
00332 public:
00333 MittelmannBndryCntrlNeum2()
00334 {}
00335
00336 virtual ~MittelmannBndryCntrlNeum2()
00337 {}
00338
00339 virtual bool InitializeProblem(Index N)
00340 {
00341 if (N<1) {
00342 printf("N has to be at least 1.");
00343 return false;
00344 }
00345 Number alpha = 0.;
00346 Number lb_y = -1e20;
00347 Number ub_y = 2.835;
00348 Number lb_u = 6.;
00349 Number ub_u = 9.;
00350 Number u_init = (ub_u+lb_u)/2.;
00351
00352 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00353 return true;
00354 }
00355 protected:
00357 virtual Number y_d_cont(Number x1, Number x2) const
00358 {
00359 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00360 }
00362 virtual Number d_cont(Number x1, Number x2, Number y) const
00363 {
00364 return 0.;
00365 }
00367 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00368 {
00369 return 0.;
00370 }
00372 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00373 {
00374 return 0.;
00375 }
00378 virtual bool d_cont_dydy_alwayszero() const
00379 {
00380 return true;
00381 }
00383 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00384 {
00385 return u - y*y;
00386 }
00388 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00389 {
00390 return - 2.*y;
00391 }
00393 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00394 {
00395 return 1.;
00396 }
00398 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00399 {
00400 return -2.;
00401 }
00404 virtual bool b_cont_dydy_alwayszero() const
00405 {
00406 return false;
00407 }
00408 private:
00411 MittelmannBndryCntrlNeum2(const MittelmannBndryCntrlNeum2&);
00412 MittelmannBndryCntrlNeum2& operator=(const MittelmannBndryCntrlNeum2&);
00414 };
00415
00417 class MittelmannBndryCntrlNeum3 : public MittelmannBndryCntrlNeumBase
00418 {
00419 public:
00420 MittelmannBndryCntrlNeum3()
00421 {}
00422
00423 virtual ~MittelmannBndryCntrlNeum3()
00424 {}
00425
00426 virtual bool InitializeProblem(Index N)
00427 {
00428 if (N<1) {
00429 printf("N has to be at least 1.");
00430 return false;
00431 }
00432 Number alpha = 0.01;
00433 Number lb_y = -1e20;
00434 Number ub_y = 2.7;
00435 Number lb_u = 1.8;
00436 Number ub_u = 2.5;
00437 Number u_init = (ub_u+lb_u)/2.;
00438
00439 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00440 return true;
00441 }
00442 protected:
00444 virtual Number y_d_cont(Number x1, Number x2) const
00445 {
00446 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00447 }
00449 virtual Number d_cont(Number x1, Number x2, Number y) const
00450 {
00451 return y*y*y-y;
00452 }
00454 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00455 {
00456 return 3.*y*y-1.;
00457 }
00459 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00460 {
00461 return 6.*y;
00462 }
00465 virtual bool d_cont_dydy_alwayszero() const
00466 {
00467 return false;
00468 }
00470 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00471 {
00472 return u;
00473 }
00475 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00476 {
00477 return 0.;
00478 }
00480 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00481 {
00482 return 1.;
00483 }
00485 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00486 {
00487 return 0.;
00488 }
00491 virtual bool b_cont_dydy_alwayszero() const
00492 {
00493 return true;
00494 }
00495 private:
00498 MittelmannBndryCntrlNeum3(const MittelmannBndryCntrlNeum3&);
00499 MittelmannBndryCntrlNeum3& operator=(const MittelmannBndryCntrlNeum3&);
00501 };
00502
00504 class MittelmannBndryCntrlNeum4 : public MittelmannBndryCntrlNeumBase
00505 {
00506 public:
00507 MittelmannBndryCntrlNeum4()
00508 {}
00509
00510 virtual ~MittelmannBndryCntrlNeum4()
00511 {}
00512
00513 virtual bool InitializeProblem(Index N)
00514 {
00515 if (N<1) {
00516 printf("N has to be at least 1.");
00517 return false;
00518 }
00519 Number alpha = 0.;
00520 Number lb_y = -1e20;
00521 Number ub_y = 2.7;
00522 Number lb_u = 1.8;
00523 Number ub_u = 2.5;
00524 Number u_init = (ub_u+lb_u)/2.;
00525
00526 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00527 return true;
00528 }
00529 protected:
00531 virtual Number y_d_cont(Number x1, Number x2) const
00532 {
00533 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00534 }
00536 virtual Number d_cont(Number x1, Number x2, Number y) const
00537 {
00538 return y*y*y-y;
00539 }
00541 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00542 {
00543 return 3.*y*y-1.;
00544 }
00546 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00547 {
00548 return 6.*y;
00549 }
00552 virtual bool d_cont_dydy_alwayszero() const
00553 {
00554 return false;
00555 }
00557 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00558 {
00559 return u;
00560 }
00562 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00563 {
00564 return 0.;
00565 }
00567 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00568 {
00569 return 1.;
00570 }
00572 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00573 {
00574 return 0.;
00575 }
00578 virtual bool b_cont_dydy_alwayszero() const
00579 {
00580 return true;
00581 }
00582 private:
00585 MittelmannBndryCntrlNeum4(const MittelmannBndryCntrlNeum4&);
00586 MittelmannBndryCntrlNeum4& operator=(const MittelmannBndryCntrlNeum4&);
00588 };
00589
00590 #endif