00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #ifndef __MITTELMANNDISTRCNTRLNEUMA_HPP__
00011 #define __MITTELMANNDISTRCNTRLNEUMA_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
00034 class MittelmannDistCntrlNeumABase : public RegisteredTNLP
00035 {
00036 public:
00039 MittelmannDistCntrlNeumABase();
00040
00042 virtual ~MittelmannDistCntrlNeumABase();
00043
00047 virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00048 Index& nnz_h_lag, IndexStyleEnum& index_style);
00049
00051 virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00052 Index m, Number* g_l, Number* g_u);
00053
00055 virtual bool get_starting_point(Index n, bool init_x, Number* x,
00056 bool init_z, Number* z_L, Number* z_U,
00057 Index m, bool init_lambda,
00058 Number* lambda);
00059
00061 virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00062
00064 virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00065
00067 virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00068
00073 virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00074 Index m, Index nele_jac, Index* iRow, Index *jCol,
00075 Number* values);
00076
00081 virtual bool eval_h(Index n, const Number* x, bool new_x,
00082 Number obj_factor, Index m, const Number* lambda,
00083 bool new_lambda, Index nele_hess, Index* iRow,
00084 Index* jCol, Number* values);
00085
00087
00089 virtual bool get_scaling_parameters(Number& obj_scaling,
00090 bool& use_x_scaling, Index n,
00091 Number* x_scaling,
00092 bool& use_g_scaling, Index m,
00093 Number* g_scaling);
00094
00099 virtual void finalize_solution(SolverReturn status,
00100 Index n, const Number* x, const Number* z_L, const Number* z_U,
00101 Index m, const Number* g, const Number* lambda,
00102 Number obj_value);
00104
00105 protected:
00109 void SetBaseParameters(Index N, Number lb_y,
00110 Number ub_y, Number lb_u, Number ub_u,
00111 Number b_0j, Number b_1j, Number b_i0, Number b_i1,
00112 Number u_init);
00113
00117 virtual Number y_d_cont(Number x1, Number x2) const =0;
00119 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const =0;
00121 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00123 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00125 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00128 virtual bool fint_cont_dydy_alwayszero() const =0;
00130 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00133 virtual bool fint_cont_dudu_alwayszero() const =0;
00135 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00138 virtual bool fint_cont_dydu_alwayszero() const =0;
00140 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const =0;
00142 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00144 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00146 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00149 virtual bool d_cont_dydy_alwayszero() const =0;
00151 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00154 virtual bool d_cont_dudu_alwayszero() const =0;
00156 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00159 virtual bool d_cont_dydu_alwayszero() const =0;
00161
00162 private:
00174 MittelmannDistCntrlNeumABase(const MittelmannDistCntrlNeumABase&);
00175 MittelmannDistCntrlNeumABase& operator=(const MittelmannDistCntrlNeumABase&);
00177
00181 Index N_;
00183 Number h_;
00185 Number hh_;
00187 Number lb_y_;
00189 Number ub_y_;
00191 Number lb_u_;
00193 Number ub_u_;
00196 Number b_0j_;
00199 Number b_1j_;
00202 Number b_i0_;
00205 Number b_i1_;
00207 Number u_init_;
00209 Number* y_d_;
00211
00216 inline Index y_index(Index i, Index j) const
00217 {
00218 return j + (N_+2)*i;
00219 }
00222 inline Index u_index(Index i, Index j) const
00223 {
00224 return (N_+2)*(N_+2) + (j-1) + (N_)*(i-1);
00225 }
00228 inline Index pde_index(Index i, Index j) const
00229 {
00230 return (j-1) + N_*(i-1);
00231 }
00233 inline Number x1_grid(Index i) const
00234 {
00235 return h_*(Number)i;
00236 }
00238 inline Number x2_grid(Index i) const
00239 {
00240 return h_*(Number)i;
00241 }
00243 };
00244
00246 class MittelmannDistCntrlNeumA1 : public MittelmannDistCntrlNeumABase
00247 {
00248 public:
00249 MittelmannDistCntrlNeumA1()
00250 :
00251 pi_(4.*atan(1.)),
00252 alpha_(0.001)
00253 {}
00254
00255 virtual ~MittelmannDistCntrlNeumA1()
00256 {}
00257
00258 virtual bool InitializeProblem(Index N)
00259 {
00260 if (N<1) {
00261 printf("N has to be at least 1.");
00262 return false;
00263 }
00264 Number lb_y = -1e20;
00265 Number ub_y = 0.371;
00266 Number lb_u = -8.;
00267 Number ub_u = 9.;
00268 Number b_0j = 1.;
00269 Number b_1j = 1.;
00270 Number b_i0 = 1.;
00271 Number b_i1 = 1.;
00272 Number u_init = (ub_u+lb_u)/2.;
00273
00274 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00275 return true;
00276 }
00277 protected:
00279 virtual Number y_d_cont(Number x1, Number x2) const
00280 {
00281 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00282 }
00284 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00285 {
00286 Number diff_y = y-y_d_cont(x1,x2);
00287 return 0.5*(diff_y*diff_y + alpha_*u*u);
00288 }
00290 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00291 {
00292 return y-y_d_cont(x1,x2);
00293 }
00294
00296 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00297 {
00298 return alpha_*u;
00299 }
00301 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00302 {
00303 return 1.;
00304 }
00307 virtual bool fint_cont_dydy_alwayszero() const
00308 {
00309 return false;
00310 }
00312 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00313 {
00314 return alpha_;
00315 }
00318 virtual bool fint_cont_dudu_alwayszero() const
00319 {
00320 return false;
00321 }
00323 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00324 {
00325 return 0.;
00326 }
00329 virtual bool fint_cont_dydu_alwayszero() const
00330 {
00331 return true;
00332 }
00334 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00335 {
00336 return -exp(y) - u;
00337 }
00339 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00340 {
00341 return -exp(y);
00342 }
00344 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00345 {
00346 return -1.;
00347 }
00349 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00350 {
00351 return -exp(y);
00352 }
00355 virtual bool d_cont_dydy_alwayszero() const
00356 {
00357 return false;
00358 }
00360 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00361 {
00362 return 0.;
00363 }
00366 virtual bool d_cont_dudu_alwayszero() const
00367 {
00368 return true;
00369 }
00371 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00372 {
00373 return 0.;
00374 }
00377 virtual bool d_cont_dydu_alwayszero() const
00378 {
00379 return true;
00380 }
00381 private:
00384 MittelmannDistCntrlNeumA1(const MittelmannDistCntrlNeumA1&);
00385 MittelmannDistCntrlNeumA1& operator=(const MittelmannDistCntrlNeumA1&);
00387
00388 const Number pi_;
00390 const Number alpha_;
00391 };
00392
00394 class MittelmannDistCntrlNeumA2 : public MittelmannDistCntrlNeumABase
00395 {
00396 public:
00397 MittelmannDistCntrlNeumA2()
00398 :
00399 pi_(4.*atan(1.))
00400 {}
00401
00402 virtual ~MittelmannDistCntrlNeumA2()
00403 {}
00404
00405 virtual bool InitializeProblem(Index N)
00406 {
00407 if (N<1) {
00408 printf("N has to be at least 1.");
00409 return false;
00410 }
00411 Number lb_y = -1e20;
00412 Number ub_y = 0.371;
00413 Number lb_u = -8.;
00414 Number ub_u = 9.;
00415 Number b_0j = 1.;
00416 Number b_1j = 1.;
00417 Number b_i0 = 1.;
00418 Number b_i1 = 1.;
00419 Number u_init = (ub_u+lb_u)/2.;
00420
00421 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00422 return true;
00423 }
00424 protected:
00426 virtual Number y_d_cont(Number x1, Number x2) const
00427 {
00428 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00429 }
00431 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00432 {
00433 Number diff_y = y-y_d_cont(x1,x2);
00434 return 0.5*diff_y*diff_y;
00435 }
00437 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00438 {
00439 return y-y_d_cont(x1,x2);
00440 }
00441
00443 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00444 {
00445 return 0.;
00446 }
00448 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00449 {
00450 return 1.;
00451 }
00454 virtual bool fint_cont_dydy_alwayszero() const
00455 {
00456 return false;
00457 }
00459 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00460 {
00461 return 0.;
00462 }
00465 virtual bool fint_cont_dudu_alwayszero() const
00466 {
00467 return true;
00468 }
00470 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00471 {
00472 return 0.;
00473 }
00476 virtual bool fint_cont_dydu_alwayszero() const
00477 {
00478 return true;
00479 }
00481 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00482 {
00483 return -exp(y) - u;
00484 }
00486 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00487 {
00488 return -exp(y);
00489 }
00491 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00492 {
00493 return -1.;
00494 }
00496 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00497 {
00498 return -exp(y);
00499 }
00502 virtual bool d_cont_dydy_alwayszero() const
00503 {
00504 return false;
00505 }
00507 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00508 {
00509 return 0.;
00510 }
00513 virtual bool d_cont_dudu_alwayszero() const
00514 {
00515 return true;
00516 }
00518 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00519 {
00520 return 0.;
00521 }
00524 virtual bool d_cont_dydu_alwayszero() const
00525 {
00526 return true;
00527 }
00528 private:
00531 MittelmannDistCntrlNeumA2(const MittelmannDistCntrlNeumA2&);
00532 MittelmannDistCntrlNeumA2& operator=(const MittelmannDistCntrlNeumA2&);
00534
00535 const Number pi_;
00536 };
00537
00539 class MittelmannDistCntrlNeumA3 : public MittelmannDistCntrlNeumABase
00540 {
00541 public:
00542 MittelmannDistCntrlNeumA3()
00543 :
00544 pi_(4.*atan(1.)),
00545 M_(1.),
00546 K_(0.8),
00547 b_(1.)
00548 {}
00549
00550 virtual ~MittelmannDistCntrlNeumA3()
00551 {}
00552
00553 virtual bool InitializeProblem(Index N)
00554 {
00555 if (N<1) {
00556 printf("N has to be at least 1.");
00557 return false;
00558 }
00559 Number lb_y = 3.;
00560 Number ub_y = 6.09;
00561 Number lb_u = 1.4;
00562 Number ub_u = 1.6;
00563 Number b_0j = 1.;
00564 Number b_1j = 0.;
00565 Number b_i0 = 1.;
00566 Number b_i1 = 0.;
00567 Number u_init = (ub_u+lb_u)/2.;
00568
00569 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00570 return true;
00571 }
00572 protected:
00574 virtual Number y_d_cont(Number x1, Number x2) const
00575 {
00576 return 6.;
00577 }
00579 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00580 {
00581 return u*(M_*u - K_*y);
00582 }
00584 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00585 {
00586 return -K_*u;
00587 }
00588
00590 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00591 {
00592 return 2.*M_*u - K_*y;
00593 }
00595 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00596 {
00597 return 0.;
00598 }
00601 virtual bool fint_cont_dydy_alwayszero() const
00602 {
00603 return true;
00604 }
00606 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00607 {
00608 return 2.*M_;
00609 }
00612 virtual bool fint_cont_dudu_alwayszero() const
00613 {
00614 return false;
00615 }
00617 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00618 {
00619 return -K_;
00620 }
00623 virtual bool fint_cont_dydu_alwayszero() const
00624 {
00625 return false;
00626 }
00628 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00629 {
00630 return y*(u + b_*y - a(x1,x2));
00631 }
00633 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00634 {
00635 return (u + 2.*b_*y -a(x1,x2));
00636 }
00638 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00639 {
00640 return y;
00641 }
00643 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00644 {
00645 return 2.*b_;
00646 }
00649 virtual bool d_cont_dydy_alwayszero() const
00650 {
00651 return false;
00652 }
00654 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00655 {
00656 return 0.;
00657 }
00660 virtual bool d_cont_dudu_alwayszero() const
00661 {
00662 return true;
00663 }
00665 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00666 {
00667 return 1.;
00668 }
00671 virtual bool d_cont_dydu_alwayszero() const
00672 {
00673 return false;
00674 }
00675 private:
00678 MittelmannDistCntrlNeumA3(const MittelmannDistCntrlNeumA3&);
00679 MittelmannDistCntrlNeumA3& operator=(const MittelmannDistCntrlNeumA3&);
00681
00682 const Number pi_;
00683
00685 const Number M_;
00686 const Number K_;
00687 const Number b_;
00689
00690 inline Number a(Number x1, Number x2) const
00691 {
00692 return 7. + 4.*sin(2.*pi_*x1*x2);
00693 }
00694 };
00695
00696 #endif