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_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
00050 class MittelmannDistCntrlNeumABase : public RegisteredTNLP
00051 {
00052 public:
00055 MittelmannDistCntrlNeumABase();
00056
00058 virtual ~MittelmannDistCntrlNeumABase();
00059
00063 virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00064 Index& nnz_h_lag, IndexStyleEnum& index_style);
00065
00067 virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00068 Index m, Number* g_l, Number* g_u);
00069
00071 virtual bool get_starting_point(Index n, bool init_x, Number* x,
00072 bool init_z, Number* z_L, Number* z_U,
00073 Index m, bool init_lambda,
00074 Number* lambda);
00075
00077 virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00078
00080 virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00081
00083 virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00084
00089 virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00090 Index m, Index nele_jac, Index* iRow, Index *jCol,
00091 Number* values);
00092
00097 virtual bool eval_h(Index n, const Number* x, bool new_x,
00098 Number obj_factor, Index m, const Number* lambda,
00099 bool new_lambda, Index nele_hess, Index* iRow,
00100 Index* jCol, Number* values);
00101
00103
00105 virtual bool get_scaling_parameters(Number& obj_scaling,
00106 bool& use_x_scaling, Index n,
00107 Number* x_scaling,
00108 bool& use_g_scaling, Index m,
00109 Number* g_scaling);
00110
00115 virtual void finalize_solution(SolverReturn status,
00116 Index n, const Number* x, const Number* z_L, const Number* z_U,
00117 Index m, const Number* g, const Number* lambda,
00118 Number obj_value,
00119 const IpoptData* ip_data,
00120 IpoptCalculatedQuantities* ip_cq);
00122
00123 protected:
00127 void SetBaseParameters(Index N, Number lb_y,
00128 Number ub_y, Number lb_u, Number ub_u,
00129 Number b_0j, Number b_1j, Number b_i0, Number b_i1,
00130 Number u_init);
00131
00135 virtual Number y_d_cont(Number x1, Number x2) const =0;
00137 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const =0;
00139 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00141 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00143 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00146 virtual bool fint_cont_dydy_alwayszero() const =0;
00148 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00151 virtual bool fint_cont_dudu_alwayszero() const =0;
00153 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00156 virtual bool fint_cont_dydu_alwayszero() const =0;
00158 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const =0;
00160 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const =0;
00162 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const =0;
00164 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const =0;
00167 virtual bool d_cont_dydy_alwayszero() const =0;
00169 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const =0;
00172 virtual bool d_cont_dudu_alwayszero() const =0;
00174 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const =0;
00177 virtual bool d_cont_dydu_alwayszero() const =0;
00179
00180 private:
00192 MittelmannDistCntrlNeumABase(const MittelmannDistCntrlNeumABase&);
00193 MittelmannDistCntrlNeumABase& operator=(const MittelmannDistCntrlNeumABase&);
00195
00199 Index N_;
00201 Number h_;
00203 Number hh_;
00205 Number lb_y_;
00207 Number ub_y_;
00209 Number lb_u_;
00211 Number ub_u_;
00214 Number b_0j_;
00217 Number b_1j_;
00220 Number b_i0_;
00223 Number b_i1_;
00225 Number u_init_;
00227 Number* y_d_;
00229
00234 inline Index y_index(Index i, Index j) const
00235 {
00236 return j + (N_+2)*i;
00237 }
00240 inline Index u_index(Index i, Index j) const
00241 {
00242 return (N_+2)*(N_+2) + (j-1) + (N_)*(i-1);
00243 }
00246 inline Index pde_index(Index i, Index j) const
00247 {
00248 return (j-1) + N_*(i-1);
00249 }
00251 inline Number x1_grid(Index i) const
00252 {
00253 return h_*(Number)i;
00254 }
00256 inline Number x2_grid(Index i) const
00257 {
00258 return h_*(Number)i;
00259 }
00261 };
00262
00264 class MittelmannDistCntrlNeumA1 : public MittelmannDistCntrlNeumABase
00265 {
00266 public:
00267 MittelmannDistCntrlNeumA1()
00268 :
00269 pi_(4.*atan(1.)),
00270 alpha_(0.001)
00271 {}
00272
00273 virtual ~MittelmannDistCntrlNeumA1()
00274 {}
00275
00276 virtual bool InitializeProblem(Index N)
00277 {
00278 if (N<1) {
00279 printf("N has to be at least 1.");
00280 return false;
00281 }
00282 Number lb_y = -1e20;
00283 Number ub_y = 0.371;
00284 Number lb_u = -8.;
00285 Number ub_u = 9.;
00286 Number b_0j = 1.;
00287 Number b_1j = 1.;
00288 Number b_i0 = 1.;
00289 Number b_i1 = 1.;
00290 Number u_init = (ub_u+lb_u)/2.;
00291
00292 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00293 return true;
00294 }
00295 protected:
00297 virtual Number y_d_cont(Number x1, Number x2) const
00298 {
00299 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00300 }
00302 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00303 {
00304 Number diff_y = y-y_d_cont(x1,x2);
00305 return 0.5*(diff_y*diff_y + alpha_*u*u);
00306 }
00308 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00309 {
00310 return y-y_d_cont(x1,x2);
00311 }
00312
00314 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00315 {
00316 return alpha_*u;
00317 }
00319 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00320 {
00321 return 1.;
00322 }
00325 virtual bool fint_cont_dydy_alwayszero() const
00326 {
00327 return false;
00328 }
00330 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00331 {
00332 return alpha_;
00333 }
00336 virtual bool fint_cont_dudu_alwayszero() const
00337 {
00338 return false;
00339 }
00341 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00342 {
00343 return 0.;
00344 }
00347 virtual bool fint_cont_dydu_alwayszero() const
00348 {
00349 return true;
00350 }
00352 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00353 {
00354 return -exp(y) - u;
00355 }
00357 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00358 {
00359 return -exp(y);
00360 }
00362 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00363 {
00364 return -1.;
00365 }
00367 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00368 {
00369 return -exp(y);
00370 }
00373 virtual bool d_cont_dydy_alwayszero() const
00374 {
00375 return false;
00376 }
00378 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00379 {
00380 return 0.;
00381 }
00384 virtual bool d_cont_dudu_alwayszero() const
00385 {
00386 return true;
00387 }
00389 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00390 {
00391 return 0.;
00392 }
00395 virtual bool d_cont_dydu_alwayszero() const
00396 {
00397 return true;
00398 }
00399 private:
00402 MittelmannDistCntrlNeumA1(const MittelmannDistCntrlNeumA1&);
00403 MittelmannDistCntrlNeumA1& operator=(const MittelmannDistCntrlNeumA1&);
00405
00406 const Number pi_;
00408 const Number alpha_;
00409 };
00410
00412 class MittelmannDistCntrlNeumA2 : public MittelmannDistCntrlNeumABase
00413 {
00414 public:
00415 MittelmannDistCntrlNeumA2()
00416 :
00417 pi_(4.*atan(1.))
00418 {}
00419
00420 virtual ~MittelmannDistCntrlNeumA2()
00421 {}
00422
00423 virtual bool InitializeProblem(Index N)
00424 {
00425 if (N<1) {
00426 printf("N has to be at least 1.");
00427 return false;
00428 }
00429 Number lb_y = -1e20;
00430 Number ub_y = 0.371;
00431 Number lb_u = -8.;
00432 Number ub_u = 9.;
00433 Number b_0j = 1.;
00434 Number b_1j = 1.;
00435 Number b_i0 = 1.;
00436 Number b_i1 = 1.;
00437 Number u_init = (ub_u+lb_u)/2.;
00438
00439 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00440 return true;
00441 }
00442 protected:
00444 virtual Number y_d_cont(Number x1, Number x2) const
00445 {
00446 return sin(2.*pi_*x1)*sin(2.*pi_*x2);
00447 }
00449 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00450 {
00451 Number diff_y = y-y_d_cont(x1,x2);
00452 return 0.5*diff_y*diff_y;
00453 }
00455 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00456 {
00457 return y-y_d_cont(x1,x2);
00458 }
00459
00461 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00462 {
00463 return 0.;
00464 }
00466 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00467 {
00468 return 1.;
00469 }
00472 virtual bool fint_cont_dydy_alwayszero() const
00473 {
00474 return false;
00475 }
00477 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00478 {
00479 return 0.;
00480 }
00483 virtual bool fint_cont_dudu_alwayszero() const
00484 {
00485 return true;
00486 }
00488 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00489 {
00490 return 0.;
00491 }
00494 virtual bool fint_cont_dydu_alwayszero() const
00495 {
00496 return true;
00497 }
00499 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00500 {
00501 return -exp(y) - u;
00502 }
00504 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00505 {
00506 return -exp(y);
00507 }
00509 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00510 {
00511 return -1.;
00512 }
00514 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00515 {
00516 return -exp(y);
00517 }
00520 virtual bool d_cont_dydy_alwayszero() const
00521 {
00522 return false;
00523 }
00525 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00526 {
00527 return 0.;
00528 }
00531 virtual bool d_cont_dudu_alwayszero() const
00532 {
00533 return true;
00534 }
00536 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00537 {
00538 return 0.;
00539 }
00542 virtual bool d_cont_dydu_alwayszero() const
00543 {
00544 return true;
00545 }
00546 private:
00549 MittelmannDistCntrlNeumA2(const MittelmannDistCntrlNeumA2&);
00550 MittelmannDistCntrlNeumA2& operator=(const MittelmannDistCntrlNeumA2&);
00552
00553 const Number pi_;
00554 };
00555
00557 class MittelmannDistCntrlNeumA3 : public MittelmannDistCntrlNeumABase
00558 {
00559 public:
00560 MittelmannDistCntrlNeumA3()
00561 :
00562 pi_(4.*atan(1.)),
00563 M_(1.),
00564 K_(0.8),
00565 b_(1.)
00566 {}
00567
00568 virtual ~MittelmannDistCntrlNeumA3()
00569 {}
00570
00571 virtual bool InitializeProblem(Index N)
00572 {
00573 if (N<1) {
00574 printf("N has to be at least 1.");
00575 return false;
00576 }
00577 Number lb_y = 3.;
00578 Number ub_y = 6.09;
00579 Number lb_u = 1.4;
00580 Number ub_u = 1.6;
00581 Number b_0j = 1.;
00582 Number b_1j = 0.;
00583 Number b_i0 = 1.;
00584 Number b_i1 = 0.;
00585 Number u_init = (ub_u+lb_u)/2.;
00586
00587 SetBaseParameters(N, lb_y, ub_y, lb_u, ub_u, b_0j, b_1j, b_i0, b_i1, u_init);
00588 return true;
00589 }
00590 protected:
00592 virtual Number y_d_cont(Number x1, Number x2) const
00593 {
00594 return 6.;
00595 }
00597 virtual Number fint_cont(Number x1, Number x2, Number y, Number u) const
00598 {
00599 return u*(M_*u - K_*y);
00600 }
00602 virtual Number fint_cont_dy(Number x1, Number x2, Number y, Number u) const
00603 {
00604 return -K_*u;
00605 }
00606
00608 virtual Number fint_cont_du(Number x1, Number x2, Number y, Number u) const
00609 {
00610 return 2.*M_*u - K_*y;
00611 }
00613 virtual Number fint_cont_dydy(Number x1, Number x2, Number y, Number u) const
00614 {
00615 return 0.;
00616 }
00619 virtual bool fint_cont_dydy_alwayszero() const
00620 {
00621 return true;
00622 }
00624 virtual Number fint_cont_dudu(Number x1, Number x2, Number y, Number u) const
00625 {
00626 return 2.*M_;
00627 }
00630 virtual bool fint_cont_dudu_alwayszero() const
00631 {
00632 return false;
00633 }
00635 virtual Number fint_cont_dydu(Number x1, Number x2, Number y, Number u) const
00636 {
00637 return -K_;
00638 }
00641 virtual bool fint_cont_dydu_alwayszero() const
00642 {
00643 return false;
00644 }
00646 virtual Number d_cont(Number x1, Number x2, Number y, Number u) const
00647 {
00648 return y*(u + b_*y - a(x1,x2));
00649 }
00651 virtual Number d_cont_dy(Number x1, Number x2, Number y, Number u) const
00652 {
00653 return (u + 2.*b_*y -a(x1,x2));
00654 }
00656 virtual Number d_cont_du(Number x1, Number x2, Number y, Number u) const
00657 {
00658 return y;
00659 }
00661 virtual Number d_cont_dydy(Number x1, Number x2, Number y, Number u) const
00662 {
00663 return 2.*b_;
00664 }
00667 virtual bool d_cont_dydy_alwayszero() const
00668 {
00669 return false;
00670 }
00672 virtual Number d_cont_dudu(Number x1, Number x2, Number y, Number u) const
00673 {
00674 return 0.;
00675 }
00678 virtual bool d_cont_dudu_alwayszero() const
00679 {
00680 return true;
00681 }
00683 virtual Number d_cont_dydu(Number x1, Number x2, Number y, Number u) const
00684 {
00685 return 1.;
00686 }
00689 virtual bool d_cont_dydu_alwayszero() const
00690 {
00691 return false;
00692 }
00693 private:
00696 MittelmannDistCntrlNeumA3(const MittelmannDistCntrlNeumA3&);
00697 MittelmannDistCntrlNeumA3& operator=(const MittelmannDistCntrlNeumA3&);
00699
00700 const Number pi_;
00701
00703 const Number M_;
00704 const Number K_;
00705 const Number b_;
00707
00708 inline Number a(Number x1, Number x2) const
00709 {
00710 return 7. + 4.*sin(2.*pi_*x1*x2);
00711 }
00712 };
00713
00714 #endif