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,
00103 const IpoptData* ip_data,
00104 IpoptCalculatedQuantities* ip_cq);
00106
00107 protected:
00111 void SetBaseParameters(Index N, Number lb_y,
00112 Number ub_y, Number lb_u, Number ub_u,
00113 Number b_0j, Number b_1j, Number b_i0, Number b_i1,
00114 Number u_init);
00115
00119 virtual Number y_d_cont(Number x1, Number x2) const =0;
00121 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const =0;
00123 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00125 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00127 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00130 virtual bool fint_cont_dydy_alwayszero() const =0;
00132 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00135 virtual bool fint_cont_dudu_alwayszero() const =0;
00137 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00140 virtual bool fint_cont_dydu_alwayszero() const =0;
00142 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const =0;
00144 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00146 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00148 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00151 virtual bool d_cont_dydy_alwayszero() const =0;
00153 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00156 virtual bool d_cont_dudu_alwayszero() const =0;
00158 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00161 virtual bool d_cont_dydu_alwayszero() const =0;
00163
00164 private:
00176 MittelmannDistCntrlNeumABase(const MittelmannDistCntrlNeumABase&);
00177 MittelmannDistCntrlNeumABase& operator=(const MittelmannDistCntrlNeumABase&);
00179
00183 Index N_;
00185 Number h_;
00187 Number hh_;
00189 Number lb_y_;
00191 Number ub_y_;
00193 Number lb_u_;
00195 Number ub_u_;
00198 Number b_0j_;
00201 Number b_1j_;
00204 Number b_i0_;
00207 Number b_i1_;
00209 Number u_init_;
00211 Number* y_d_;
00213
00218 inline Index y_index(Index i, Index j) const
00219 {
00220 return j + (N_+2)*i;
00221 }
00224 inline Index u_index(Index i, Index j) const
00225 {
00226 return (N_+2)*(N_+2) + (j-1) + (N_)*(i-1);
00227 }
00230 inline Index pde_index(Index i, Index j) const
00231 {
00232 return (j-1) + N_*(i-1);
00233 }
00235 inline Number x1_grid(Index i) const
00236 {
00237 return h_*(Number)i;
00238 }
00240 inline Number x2_grid(Index i) const
00241 {
00242 return h_*(Number)i;
00243 }
00245 };
00246
00248 class MittelmannDistCntrlNeumA1 : public MittelmannDistCntrlNeumABase
00249 {
00250 public:
00251 MittelmannDistCntrlNeumA1()
00252 :
00253 pi_(4.*atan(1.)),
00254 alpha_(0.001)
00255 {}
00256
00257 virtual ~MittelmannDistCntrlNeumA1()
00258 {}
00259
00260 virtual bool InitializeProblem(Index N)
00261 {
00262 if (N<1) {
00263 printf("N has to be at least 1.");
00264 return false;
00265 }
00266 Number lb_y = -1e20;
00267 Number ub_y = 0.371;
00268 Number lb_u = -8.;
00269 Number ub_u = 9.;
00270 Number b_0j = 1.;
00271 Number b_1j = 1.;
00272 Number b_i0 = 1.;
00273 Number b_i1 = 1.;
00274 Number u_init = (ub_u+lb_u)/2.;
00275
00276 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00277 return true;
00278 }
00279 protected:
00281 virtual Number y_d_cont(Number x1, Number x2) const
00282 {
00283 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00284 }
00286 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00287 {
00288 Number diff_y = y-y_d_cont(x1,x2);
00289 return 0.5*(diff_y*diff_y + alpha_*u*u);
00290 }
00292 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00293 {
00294 return y-y_d_cont(x1,x2);
00295 }
00296
00298 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00299 {
00300 return alpha_*u;
00301 }
00303 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00304 {
00305 return 1.;
00306 }
00309 virtual bool fint_cont_dydy_alwayszero() const
00310 {
00311 return false;
00312 }
00314 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00315 {
00316 return alpha_;
00317 }
00320 virtual bool fint_cont_dudu_alwayszero() const
00321 {
00322 return false;
00323 }
00325 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00326 {
00327 return 0.;
00328 }
00331 virtual bool fint_cont_dydu_alwayszero() const
00332 {
00333 return true;
00334 }
00336 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00337 {
00338 return -exp(y) - u;
00339 }
00341 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00342 {
00343 return -exp(y);
00344 }
00346 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00347 {
00348 return -1.;
00349 }
00351 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00352 {
00353 return -exp(y);
00354 }
00357 virtual bool d_cont_dydy_alwayszero() const
00358 {
00359 return false;
00360 }
00362 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00363 {
00364 return 0.;
00365 }
00368 virtual bool d_cont_dudu_alwayszero() const
00369 {
00370 return true;
00371 }
00373 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00374 {
00375 return 0.;
00376 }
00379 virtual bool d_cont_dydu_alwayszero() const
00380 {
00381 return true;
00382 }
00383 private:
00386 MittelmannDistCntrlNeumA1(const MittelmannDistCntrlNeumA1&);
00387 MittelmannDistCntrlNeumA1& operator=(const MittelmannDistCntrlNeumA1&);
00389
00390 const Number pi_;
00392 const Number alpha_;
00393 };
00394
00396 class MittelmannDistCntrlNeumA2 : public MittelmannDistCntrlNeumABase
00397 {
00398 public:
00399 MittelmannDistCntrlNeumA2()
00400 :
00401 pi_(4.*atan(1.))
00402 {}
00403
00404 virtual ~MittelmannDistCntrlNeumA2()
00405 {}
00406
00407 virtual bool InitializeProblem(Index N)
00408 {
00409 if (N<1) {
00410 printf("N has to be at least 1.");
00411 return false;
00412 }
00413 Number lb_y = -1e20;
00414 Number ub_y = 0.371;
00415 Number lb_u = -8.;
00416 Number ub_u = 9.;
00417 Number b_0j = 1.;
00418 Number b_1j = 1.;
00419 Number b_i0 = 1.;
00420 Number b_i1 = 1.;
00421 Number u_init = (ub_u+lb_u)/2.;
00422
00423 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00424 return true;
00425 }
00426 protected:
00428 virtual Number y_d_cont(Number x1, Number x2) const
00429 {
00430 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00431 }
00433 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00434 {
00435 Number diff_y = y-y_d_cont(x1,x2);
00436 return 0.5*diff_y*diff_y;
00437 }
00439 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00440 {
00441 return y-y_d_cont(x1,x2);
00442 }
00443
00445 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00446 {
00447 return 0.;
00448 }
00450 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00451 {
00452 return 1.;
00453 }
00456 virtual bool fint_cont_dydy_alwayszero() const
00457 {
00458 return false;
00459 }
00461 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00462 {
00463 return 0.;
00464 }
00467 virtual bool fint_cont_dudu_alwayszero() const
00468 {
00469 return true;
00470 }
00472 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00473 {
00474 return 0.;
00475 }
00478 virtual bool fint_cont_dydu_alwayszero() const
00479 {
00480 return true;
00481 }
00483 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00484 {
00485 return -exp(y) - u;
00486 }
00488 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00489 {
00490 return -exp(y);
00491 }
00493 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00494 {
00495 return -1.;
00496 }
00498 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00499 {
00500 return -exp(y);
00501 }
00504 virtual bool d_cont_dydy_alwayszero() const
00505 {
00506 return false;
00507 }
00509 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00510 {
00511 return 0.;
00512 }
00515 virtual bool d_cont_dudu_alwayszero() const
00516 {
00517 return true;
00518 }
00520 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00521 {
00522 return 0.;
00523 }
00526 virtual bool d_cont_dydu_alwayszero() const
00527 {
00528 return true;
00529 }
00530 private:
00533 MittelmannDistCntrlNeumA2(const MittelmannDistCntrlNeumA2&);
00534 MittelmannDistCntrlNeumA2& operator=(const MittelmannDistCntrlNeumA2&);
00536
00537 const Number pi_;
00538 };
00539
00541 class MittelmannDistCntrlNeumA3 : public MittelmannDistCntrlNeumABase
00542 {
00543 public:
00544 MittelmannDistCntrlNeumA3()
00545 :
00546 pi_(4.*atan(1.)),
00547 M_(1.),
00548 K_(0.8),
00549 b_(1.)
00550 {}
00551
00552 virtual ~MittelmannDistCntrlNeumA3()
00553 {}
00554
00555 virtual bool InitializeProblem(Index N)
00556 {
00557 if (N<1) {
00558 printf("N has to be at least 1.");
00559 return false;
00560 }
00561 Number lb_y = 3.;
00562 Number ub_y = 6.09;
00563 Number lb_u = 1.4;
00564 Number ub_u = 1.6;
00565 Number b_0j = 1.;
00566 Number b_1j = 0.;
00567 Number b_i0 = 1.;
00568 Number b_i1 = 0.;
00569 Number u_init = (ub_u+lb_u)/2.;
00570
00571 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00572 return true;
00573 }
00574 protected:
00576 virtual Number y_d_cont(Number x1, Number x2) const
00577 {
00578 return 6.;
00579 }
00581 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00582 {
00583 return u*(M_*u - K_*y);
00584 }
00586 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00587 {
00588 return -K_*u;
00589 }
00590
00592 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00593 {
00594 return 2.*M_*u - K_*y;
00595 }
00597 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00598 {
00599 return 0.;
00600 }
00603 virtual bool fint_cont_dydy_alwayszero() const
00604 {
00605 return true;
00606 }
00608 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00609 {
00610 return 2.*M_;
00611 }
00614 virtual bool fint_cont_dudu_alwayszero() const
00615 {
00616 return false;
00617 }
00619 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00620 {
00621 return -K_;
00622 }
00625 virtual bool fint_cont_dydu_alwayszero() const
00626 {
00627 return false;
00628 }
00630 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00631 {
00632 return y*(u + b_*y - a(x1,x2));
00633 }
00635 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00636 {
00637 return (u + 2.*b_*y -a(x1,x2));
00638 }
00640 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00641 {
00642 return y;
00643 }
00645 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00646 {
00647 return 2.*b_;
00648 }
00651 virtual bool d_cont_dydy_alwayszero() const
00652 {
00653 return false;
00654 }
00656 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00657 {
00658 return 0.;
00659 }
00662 virtual bool d_cont_dudu_alwayszero() const
00663 {
00664 return true;
00665 }
00667 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00668 {
00669 return 1.;
00670 }
00673 virtual bool d_cont_dydu_alwayszero() const
00674 {
00675 return false;
00676 }
00677 private:
00680 MittelmannDistCntrlNeumA3(const MittelmannDistCntrlNeumA3&);
00681 MittelmannDistCntrlNeumA3& operator=(const MittelmannDistCntrlNeumA3&);
00683
00684 const Number pi_;
00685
00687 const Number M_;
00688 const Number K_;
00689 const Number b_;
00691
00692 inline Number a(Number x1, Number x2) const
00693 {
00694 return 7. + 4.*sin(2.*pi_*x1*x2);
00695 }
00696 };
00697
00698 #endif