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);
00103
00104 protected:
00108 void SetBaseParameters(Index N, Number alpha, Number lb_y,
00109 Number ub_y, Number lb_u, Number ub_u,
00110 Number u_init);
00111
00115 virtual Number y_d_cont(Number x1, Number x2) const =0;
00117 virtual Number d_cont(Number x1, Number x2, Number y) const =0;
00119 virtual Number d_cont_dy(Number x1, Number x2, Number y) const =0;
00121 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const =0;
00124 virtual bool d_cont_dydy_alwayszero() const =0;
00126 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const =0;
00128 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00130 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00132 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00135 virtual bool b_cont_dydy_alwayszero() const =0;
00137
00138 private:
00150 MittelmannBndryCntrlNeumBase(const MittelmannBndryCntrlNeumBase&);
00151 MittelmannBndryCntrlNeumBase& operator=(const MittelmannBndryCntrlNeumBase&);
00153
00157 Index N_;
00159 Number h_;
00161 Number hh_;
00163 Number lb_y_;
00165 Number ub_y_;
00167 Number lb_u_;
00169 Number ub_u_;
00171 Number u_init_;
00174 Number alpha_;
00176 Number* y_d_;
00178
00183 inline Index y_index(Index i, Index j) const
00184 {
00185 return j + (N_+2)*i;
00186 }
00189 inline Index u0j_index(Index j) const
00190 {
00191 return (N_+2)*(N_+2) + j-1;
00192 }
00195 inline Index u1j_index(Index j) const
00196 {
00197 return (N_+2)*(N_+2) + N_ + j-1;
00198 }
00201 inline Index ui0_index(Index j) const
00202 {
00203 return (N_+2)*(N_+2) + 2*N_ + j-1;
00204 }
00207 inline Index ui1_index(Index j) const
00208 {
00209 return (N_+2)*(N_+2) + 3*N_ + j-1;
00210 }
00212 inline Number x1_grid(Index i) const
00213 {
00214 return h_*(Number)i;
00215 }
00217 inline Number x2_grid(Index j) const
00218 {
00219 return h_*(Number)j;
00220 }
00222 };
00223
00225 class MittelmannBndryCntrlNeum1 : public MittelmannBndryCntrlNeumBase
00226 {
00227 public:
00228 MittelmannBndryCntrlNeum1()
00229 {}
00230
00231 virtual ~MittelmannBndryCntrlNeum1()
00232 {}
00233
00234 virtual bool InitializeProblem(Index N)
00235 {
00236 if (N<1) {
00237 printf("N has to be at least 1.");
00238 return false;
00239 }
00240 Number alpha = 0.01;
00241 Number lb_y = -1e20;
00242 Number ub_y = 2.071;
00243 Number lb_u = 3.7;
00244 Number ub_u = 4.5;
00245 Number u_init = (ub_u+lb_u)/2.;
00246
00247 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00248 return true;
00249 }
00250 protected:
00252 virtual Number y_d_cont(Number x1, Number x2) const
00253 {
00254 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00255 }
00257 virtual Number d_cont(Number x1, Number x2, Number y) const
00258 {
00259 return 0.;
00260 }
00262 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00263 {
00264 return 0.;
00265 }
00267 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00268 {
00269 return 0.;
00270 }
00273 virtual bool d_cont_dydy_alwayszero() const
00274 {
00275 return true;
00276 }
00278 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00279 {
00280 return u - y*y;
00281 }
00283 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00284 {
00285 return - 2.*y;
00286 }
00288 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00289 {
00290 return 1.;
00291 }
00293 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00294 {
00295 return -2.;
00296 }
00299 virtual bool b_cont_dydy_alwayszero() const
00300 {
00301 return false;
00302 }
00303 private:
00306 MittelmannBndryCntrlNeum1(const MittelmannBndryCntrlNeum1&);
00307 MittelmannBndryCntrlNeum1& operator=(const MittelmannBndryCntrlNeum1&);
00309 };
00310
00312 class MittelmannBndryCntrlNeum2 : public MittelmannBndryCntrlNeumBase
00313 {
00314 public:
00315 MittelmannBndryCntrlNeum2()
00316 {}
00317
00318 virtual ~MittelmannBndryCntrlNeum2()
00319 {}
00320
00321 virtual bool InitializeProblem(Index N)
00322 {
00323 if (N<1) {
00324 printf("N has to be at least 1.");
00325 return false;
00326 }
00327 Number alpha = 0.;
00328 Number lb_y = -1e20;
00329 Number ub_y = 2.835;
00330 Number lb_u = 6.;
00331 Number ub_u = 9.;
00332 Number u_init = (ub_u+lb_u)/2.;
00333
00334 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00335 return true;
00336 }
00337 protected:
00339 virtual Number y_d_cont(Number x1, Number x2) const
00340 {
00341 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00342 }
00344 virtual Number d_cont(Number x1, Number x2, Number y) const
00345 {
00346 return 0.;
00347 }
00349 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00350 {
00351 return 0.;
00352 }
00354 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00355 {
00356 return 0.;
00357 }
00360 virtual bool d_cont_dydy_alwayszero() const
00361 {
00362 return true;
00363 }
00365 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00366 {
00367 return u - y*y;
00368 }
00370 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00371 {
00372 return - 2.*y;
00373 }
00375 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00376 {
00377 return 1.;
00378 }
00380 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00381 {
00382 return -2.;
00383 }
00386 virtual bool b_cont_dydy_alwayszero() const
00387 {
00388 return false;
00389 }
00390 private:
00393 MittelmannBndryCntrlNeum2(const MittelmannBndryCntrlNeum2&);
00394 MittelmannBndryCntrlNeum2& operator=(const MittelmannBndryCntrlNeum2&);
00396 };
00397
00399 class MittelmannBndryCntrlNeum3 : public MittelmannBndryCntrlNeumBase
00400 {
00401 public:
00402 MittelmannBndryCntrlNeum3()
00403 {}
00404
00405 virtual ~MittelmannBndryCntrlNeum3()
00406 {}
00407
00408 virtual bool InitializeProblem(Index N)
00409 {
00410 if (N<1) {
00411 printf("N has to be at least 1.");
00412 return false;
00413 }
00414 Number alpha = 0.01;
00415 Number lb_y = -1e20;
00416 Number ub_y = 2.7;
00417 Number lb_u = 1.8;
00418 Number ub_u = 2.5;
00419 Number u_init = (ub_u+lb_u)/2.;
00420
00421 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00422 return true;
00423 }
00424 protected:
00426 virtual Number y_d_cont(Number x1, Number x2) const
00427 {
00428 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00429 }
00431 virtual Number d_cont(Number x1, Number x2, Number y) const
00432 {
00433 return y*y*y-y;
00434 }
00436 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00437 {
00438 return 3.*y*y-1.;
00439 }
00441 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00442 {
00443 return 6.*y;
00444 }
00447 virtual bool d_cont_dydy_alwayszero() const
00448 {
00449 return false;
00450 }
00452 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00453 {
00454 return u;
00455 }
00457 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00458 {
00459 return 0.;
00460 }
00462 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00463 {
00464 return 1.;
00465 }
00467 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00468 {
00469 return 0.;
00470 }
00473 virtual bool b_cont_dydy_alwayszero() const
00474 {
00475 return true;
00476 }
00477 private:
00480 MittelmannBndryCntrlNeum3(const MittelmannBndryCntrlNeum3&);
00481 MittelmannBndryCntrlNeum3& operator=(const MittelmannBndryCntrlNeum3&);
00483 };
00484
00486 class MittelmannBndryCntrlNeum4 : public MittelmannBndryCntrlNeumBase
00487 {
00488 public:
00489 MittelmannBndryCntrlNeum4()
00490 {}
00491
00492 virtual ~MittelmannBndryCntrlNeum4()
00493 {}
00494
00495 virtual bool InitializeProblem(Index N)
00496 {
00497 if (N<1) {
00498 printf("N has to be at least 1.");
00499 return false;
00500 }
00501 Number alpha = 0.;
00502 Number lb_y = -1e20;
00503 Number ub_y = 2.7;
00504 Number lb_u = 1.8;
00505 Number ub_u = 2.5;
00506 Number u_init = (ub_u+lb_u)/2.;
00507
00508 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, u_init);
00509 return true;
00510 }
00511 protected:
00513 virtual Number y_d_cont(Number x1, Number x2) const
00514 {
00515 return 2. - 2.*(x1*(x1-1.) + x2*(x2-1.));
00516 }
00518 virtual Number d_cont(Number x1, Number x2, Number y) const
00519 {
00520 return y*y*y-y;
00521 }
00523 virtual Number d_cont_dy(Number x1, Number x2, Number y) const
00524 {
00525 return 3.*y*y-1.;
00526 }
00528 virtual Number d_cont_dydy(Number x1, Number x2, Number y) const
00529 {
00530 return 6.*y;
00531 }
00534 virtual bool d_cont_dydy_alwayszero() const
00535 {
00536 return false;
00537 }
00539 virtual Number b_cont(Number x1, Number x2, Number y, Number u) const
00540 {
00541 return u;
00542 }
00544 virtual Number b_cont_dy(Number x1, Number x2, Number y, Number u) const
00545 {
00546 return 0.;
00547 }
00549 virtual Number b_cont_du(Number x1, Number x2, Number y, Number u) const
00550 {
00551 return 1.;
00552 }
00554 virtual Number b_cont_dydy(Number x1, Number x2, Number y, Number u) const
00555 {
00556 return 0.;
00557 }
00560 virtual bool b_cont_dydy_alwayszero() const
00561 {
00562 return true;
00563 }
00564 private:
00567 MittelmannBndryCntrlNeum4(const MittelmannBndryCntrlNeum4&);
00568 MittelmannBndryCntrlNeum4& operator=(const MittelmannBndryCntrlNeum4&);
00570 };
00571
00572 #endif