00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #ifndef __MITTELMANNBNDRYCNTRLDIRI3D_27_HPP__
00012 #define __MITTELMANNBNDRYCNTRLDIRI3D_27_HPP__
00013
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
00053 class MittelmannBndryCntrlDiriBase3D_27 : public RegisteredTNLP
00054 {
00055 public:
00057 MittelmannBndryCntrlDiriBase3D_27();
00058
00060 virtual ~MittelmannBndryCntrlDiriBase3D_27();
00061
00065 virtual bool get_nlp_info(Index& n, Index& m, Index& nnz_jac_g,
00066 Index& nnz_h_lag, IndexStyleEnum& index_style);
00067
00069 virtual bool get_bounds_info(Index n, Number* x_l, Number* x_u,
00070 Index m, Number* g_l, Number* g_u);
00071
00073 virtual bool get_starting_point(Index n, bool init_x, Number* x,
00074 bool init_z, Number* z_L, Number* z_U,
00075 Index m, bool init_lambda,
00076 Number* lambda);
00077
00079 virtual bool eval_f(Index n, const Number* x, bool new_x, Number& obj_value);
00080
00082 virtual bool eval_grad_f(Index n, const Number* x, bool new_x, Number* grad_f);
00083
00085 virtual bool eval_g(Index n, const Number* x, bool new_x, Index m, Number* g);
00086
00091 virtual bool eval_jac_g(Index n, const Number* x, bool new_x,
00092 Index m, Index nele_jac, Index* iRow, Index *jCol,
00093 Number* values);
00094
00099 virtual bool eval_h(Index n, const Number* x, bool new_x,
00100 Number obj_factor, Index m, const Number* lambda,
00101 bool new_lambda, Index nele_hess, Index* iRow,
00102 Index* jCol, Number* values);
00103
00105
00107 virtual bool get_scaling_parameters(Number& obj_scaling,
00108 bool& use_x_scaling, Index n,
00109 Number* x_scaling,
00110 bool& use_g_scaling, Index m,
00111 Number* g_scaling);
00112
00117 virtual void finalize_solution(SolverReturn status,
00118 Index n, const Number* x, const Number* z_L, const Number* z_U,
00119 Index m, const Number* g, const Number* lambda,
00120 Number obj_valu,
00121 const IpoptData* ip_data,
00122 IpoptCalculatedQuantities* ip_cq);
00124
00125 protected:
00129 void SetBaseParameters(Index N, Number alpha, Number lb_y,
00130 Number ub_y, Number lb_u, Number ub_u,
00131 Number d_const, Number B, Number C);
00132
00136 virtual Number y_d_cont(Number x1, Number x2, Number x3) const =0;
00138
00139 private:
00151 MittelmannBndryCntrlDiriBase3D_27(const MittelmannBndryCntrlDiriBase3D_27&);
00152 MittelmannBndryCntrlDiriBase3D_27& operator=(const MittelmannBndryCntrlDiriBase3D_27&);
00154
00158 Index N_;
00160 Number h_;
00162 Number hh_;
00164 Number hhh_;
00166 Number lb_y_;
00168 Number ub_y_;
00170 Number lb_u_;
00172 Number ub_u_;
00174 Number d_const_;
00177 Number alpha_;
00179 Number* y_d_;
00181
00186 inline Index y_index(Index i, Index j, Index k) const
00187 {
00188 return k + (N_+2)*j + (N_+2)*(N_+2)*i;
00189 }
00192 inline Index pde_index(Index i, Index j, Index k) const
00193 {
00194 return (k-1) + N_*(j-1) + N_*N_*(i-1);
00195 }
00197 inline Number x1_grid(Index i) const
00198 {
00199 return h_*(Number)i;
00200 }
00202 inline Number x2_grid(Index i) const
00203 {
00204 return h_*(Number)i;
00205 }
00207 inline Number x3_grid(Index i) const
00208 {
00209 return h_*(Number)i;
00210 }
00212 inline Number PenObj(Number t) const
00213 {
00214 if (B_ == 0.) {
00215 return 0.5*t*t;
00216 }
00217 else if (t > B_) {
00218 return B_*B_/2. + C_*(t - B_);
00219 }
00220 else if (t < -B_) {
00221 return B_*B_/2. + C_*(-t - B_);
00222 }
00223 else {
00224 const Number t2 = t*t;
00225 const Number t4 = t2*t2;
00226 const Number t6 = t4*t2;
00227 return PenA_*t2 + PenB_*t4 + PenC_*t6;
00228 }
00229 }
00231 inline Number PenObj_1(Number t) const
00232 {
00233 if (B_ == 0.) {
00234 return t;
00235 }
00236 else if (t > B_) {
00237 return C_;
00238 }
00239 else if (t < -B_) {
00240 return -C_;
00241 }
00242 else {
00243 const Number t2 = t*t;
00244 const Number t3 = t*t2;
00245 const Number t5 = t3*t2;
00246 return 2.*PenA_*t + 4.*PenB_*t3 + 6.*PenC_*t5;
00247 }
00248 }
00250 inline Number PenObj_2(Number t) const
00251 {
00252 if (B_ == 0.) {
00253 return 1.;
00254 }
00255 else if (t > B_) {
00256 return 0.;
00257 }
00258 else if (t < -B_) {
00259 return 0.;
00260 }
00261 else {
00262 const Number t2 = t*t;
00263 const Number t4 = t2*t2;
00264 return 2.*PenA_ + 12.*PenB_*t2 + 30.*PenC_*t4;
00265 }
00266 }
00268
00271 Number B_;
00272 Number C_;
00273 Number PenA_;
00274 Number PenB_;
00275 Number PenC_;
00277 };
00278
00280 class MittelmannBndryCntrlDiri3D_27 : public MittelmannBndryCntrlDiriBase3D_27
00281 {
00282 public:
00283 MittelmannBndryCntrlDiri3D_27()
00284 {}
00285
00286 virtual ~MittelmannBndryCntrlDiri3D_27()
00287 {}
00288
00289 virtual bool InitializeProblem(Index N)
00290 {
00291 if (N<1) {
00292 printf("N has to be at least 1.");
00293 return false;
00294 }
00295 Number alpha = 1e-2;
00296 Number lb_y = -1e20;
00297 Number ub_y = 3.5;
00298 Number lb_u = 0.;
00299 Number ub_u = 10.;
00300 Number d_const = -20.;
00301 Number B = 0.;
00302 Number C = 0.;
00303 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, d_const, B, C);
00304 return true;
00305 }
00306 protected:
00308 virtual Number y_d_cont(Number x1, Number x2, Number x3) const
00309 {
00310 return 3. + 5.*(x1*(x1-1.)*x2*(x2-1.));
00311 }
00312 private:
00315 MittelmannBndryCntrlDiri3D_27(const MittelmannBndryCntrlDiri3D_27&);
00316 MittelmannBndryCntrlDiri3D_27& operator=(const MittelmannBndryCntrlDiri3D_27&);
00318
00319 };
00320
00323 class MittelmannBndryCntrlDiri3D_27BT : public MittelmannBndryCntrlDiriBase3D_27
00324 {
00325 public:
00326 MittelmannBndryCntrlDiri3D_27BT()
00327 {}
00328
00329 virtual ~MittelmannBndryCntrlDiri3D_27BT()
00330 {}
00331
00332 virtual bool InitializeProblem(Index N)
00333 {
00334 if (N<1) {
00335 printf("N has to be at least 1.");
00336 return false;
00337 }
00338 Number alpha = 1e-2;
00339 Number lb_y = -1e20;
00340 Number ub_y = 3.5;
00341 Number lb_u = 0.;
00342 Number ub_u = 10.;
00343 Number d_const = -20.;
00344 Number B = .25;
00345 Number C = 0.01;
00346 SetBaseParameters(N, alpha, lb_y, ub_y, lb_u, ub_u, d_const, B, C);
00347 return true;
00348 }
00349 protected:
00351 virtual Number y_d_cont(Number x1, Number x2, Number x3) const
00352 {
00353 return 3. + 5.*(x1*(x1-1.)*x2*(x2-1.));
00354 }
00355 private:
00358 MittelmannBndryCntrlDiri3D_27BT(const MittelmannBndryCntrlDiri3D_27BT&);
00359 MittelmannBndryCntrlDiri3D_27BT& operator=(const MittelmannBndryCntrlDiri3D_27BT&);
00361
00362 };
00363
00364 #endif