9 #ifndef __MITTELMANNPARACNTRL_HPP__ 10 #define __MITTELMANNPARACNTRL_HPP__ 17 #include "configall_system.h" 26 # error "don't have header file for math" 30 using namespace Ipopt;
58 virtual bool get_starting_point(
Index n,
bool init_x,
Number*
x,
60 Index m,
bool init_lambda,
93 bool& use_x_scaling,
Index n,
95 bool& use_g_scaling,
Index m,
110 virtual bool InitializeProblem(
Index N);
169 return jx + (Nx_+1)*it;
173 return (Nt_+1)*(Nx_+1) + it - 1;
208 typename T::ProblemSpecs p;
211 printf(
"N has to be at least 1.");
249 typename T::ProblemSpecs p;
258 if (!p.phi_dydy_always_zero()) {
272 typename T::ProblemSpecs p;
298 x_u[iy] = x_l[iy] = p.a(
x_grid(jx));
309 g_l[Nt_*(
Nx_-1) + Nt_ + i]
310 = g_u[Nt_*(
Nx_-1) + Nt_ + i]
321 Index m,
bool init_lambda,
324 DBG_ASSERT(init_x==
true && init_z==
false && init_lambda==
false);
356 use_x_scaling =
false;
357 use_g_scaling =
false;
364 bool new_x,
Number& obj_value)
375 obj_value = 0.5*
dx_*sum;
392 obj_value +=
dt_*sum;
434 typename T::ProblemSpecs p;
474 typename T::ProblemSpecs p;
478 if (values == NULL) {
492 jCol[ijac] =
y_index(jx-1,it+1);
498 jCol[ijac] =
y_index(jx+1,it+1);
521 jCol[ijac] =
y_index(Nx_-2,it);
524 jCol[ijac] =
y_index(Nx_-1,it);
543 values[ijac++] = f2 - 2.*f;
546 values[ijac++] = -f2 - 2.*f;
553 values[ijac++] = -4.;
560 values[ijac++] = -4.*f;
561 values[ijac++] = 3.*f +
beta_ + p.phi_dy(x[
y_index(Nx_,it)]);
562 values[ijac++] = -1.;
577 typename T::ProblemSpecs p;
581 if (values == NULL) {
596 if (!p.phi_dydy_always_zero()) {
606 values[ihes++] = obj_factor*0.5*
dx_;
608 values[ihes++] = obj_factor*
dx_;
610 values[ihes++] = obj_factor*0.5*
dx_;
615 values[ihes++] = obj_factor*0.5*
alpha_*
dt_;
618 if (!p.phi_dydy_always_zero()) {
619 Index ig = (Nx_-1)*Nt_ + Nt_;
621 values[ihes++] = lambda[ig++]*p.phi_dydy(x[
y_index(Nx_,it)]);
683 return sqrt2_/2.*(exp23_-exp13_);
691 return (exp1_ + expm1_)*cos(x);
703 return sqrt2_/2.*exp13_;
708 -
Min(1.,
Max(0.,(exp(t)-exp13_)/(exp23_-exp13_)));
712 return y*pow(fabs(y),3);
716 return 4.*pow(fabs(y),3);
1015 return sqrt2_/2.*(exp23_-exp13_);
1023 return (exp1_ + expm1_)*cos(x);
1035 return sqrt2_/2.*exp13_;
1039 return exp(-4.*t)/4.
1040 -
Min(1.,
Max(0.,(exp(t)-exp13_)/(exp23_-exp13_)));
1044 return -y*sin(y/10.);
1048 return -y*cos(y/10.)/10. - sin(y/10.);
1052 return y*sin(y/10.)/100.;
virtual bool get_bounds_info(Index n, Number *x_l, Number *x_u, Index m, Number *g_l, Number *g_u)
Method to return the bounds for my problem.
Number lb_u_
overall lower bound on u
Number * x
Input: Starting point Output: Optimal solution.
Index Max(Index a, Index b)
Class for all IPOPT specific calculated quantities.
Number Number Index Number Number Index Index Index index_style
indexing style for iRow & jCol, 0 for C style, 1 for Fortran style
Number ub_u_
overall upper bound on u
virtual bool get_starting_point(Index n, bool init_x, Number *x, bool init_z, Number *z_L, Number *z_U, Index m, bool init_lambda, Number *lambda)
Method to return the starting point for the algorithm.
Number phi_dydy(Number y)
Number Number Index m
Number of constraints.
Number * a_u_
Array for weighting function a_u in objective.
Number * y_T_
Array for the target profile for y in objective.
Number * a_y_
Array for weighting function a_y in objective.
Number x_grid(Index j) const
Compute the grid coordinate for given index in x direction.
Number Number * g
Values of constraint at final point (output only - ignored if set to NULL)
Number Number Index Number Number Index Index Index Eval_F_CB Eval_G_CB Eval_Grad_F_CB Eval_Jac_G_CB Eval_H_CB eval_h
Callback function for evaluating Hessian of Lagrangian function.
Number phi_dydy(Number y)
virtual bool eval_jac_g(Index n, const Number *x, bool new_x, Index m, Index nele_jac, Index *iRow, Index *jCol, Number *values)
Method to return: 1) The structure of the jacobian (if "values" is NULL) 2) The values of the jacobia...
double Number
Type of all numbers.
MittelmannParaCntrlBase()
Constructor.
Number Number Index Number Number Index Index Index Eval_F_CB Eval_G_CB Eval_Grad_F_CB eval_grad_f
Callback function for evaluating gradient of objective function.
virtual bool eval_grad_f(Index n, const Number *x, bool new_x, Number *grad_f)
Method to return the gradient of the objective.
virtual bool get_nlp_info(Index &n, Index &m, Index &nnz_jac_g, Index &nnz_h_lag, IndexStyleEnum &index_style)
Method to return some info about the nlp.
Number dx_
Step size in x direction.
Index Nx_
Number of mesh points in x direction.
bool phi_dydy_always_zero()
bool phi_dydy_always_zero()
bool phi_dydy_always_zero()
Number Number Index Number Number Index Index Index Eval_F_CB Eval_G_CB Eval_Grad_F_CB Eval_Jac_G_CB eval_jac_g
Callback function for evaluating Jacobian of constraint functions.
SolverReturn
enum for the return from the optimize algorithm (obviously we need to add more)
Number phi_dydy(Number y)
Number Number Index Number Number Index nele_jac
Number of non-zero elements in constraint Jacobian.
Index u_index(Index it) const
virtual bool eval_h(Index n, const Number *x, bool new_x, Number obj_factor, Index m, const Number *lambda, bool new_lambda, Index nele_hess, Index *iRow, Index *jCol, Number *values)
Method to return: 1) The structure of the hessian of the lagrangian (if "values" is NULL) 2) The valu...
Index Min(Index a, Index b)
Number phi_dydy(Number y)
Class to organize all the data required by the algorithm.
Number Number Index Number Number Index Index Index Eval_F_CB Eval_G_CB eval_g
Callback function for evaluating constraint functions.
int Index
Type of all indices of vectors, matrices etc.
Index y_index(Index jx, Index it) const
Translation of mesh point indices to NLP variable indices for y(x_ij)
virtual bool get_scaling_parameters(Number &obj_scaling, bool &use_x_scaling, Index n, Number *x_scaling, bool &use_g_scaling, Index m, Number *g_scaling)
Method for returning scaling parameters.
virtual void finalize_solution(SolverReturn status, Index n, const Number *x, const Number *z_L, const Number *z_U, Index m, const Number *g, const Number *lambda, Number obj_value, const IpoptData *ip_data, IpoptCalculatedQuantities *ip_cq)
This method is called after the optimization, and could write an output file with the optimal profile...
virtual ~MittelmannParaCntrlBase()
Default destructor.
Number Number * x_scaling
bool phi_dydy_always_zero()
bool phi_dydy_always_zero()
Number phi_dydy(Number y)
Number Number Index Number Number Index Index nele_hess
Number of non-zero elements in Hessian of Lagrangian.
Class implemented the NLP discretization of.
Number t_grid(Index i) const
Compute the grid coordinate for given index in t direction.
Number l_
Upper bound on x.
Number Number Number * g_scaling
Number beta_
Weighting parameter in PDE.
IndexStyleEnum
overload this method to return the number of variables and constraints, and the number of non-zeros i...
virtual bool InitializeProblem(Index N)
Initialize internal parameters, where N is a parameter determining the problme size.
Number ub_y_
overall upper bound on y
Base class for parabolic and elliptic control problems, as formulated by Hans Mittelmann as problem (...
Number T_
Upper bound on t.
Number alpha_
Weighting parameter for the control target deviation functional in the objective. ...
virtual bool eval_g(Index n, const Number *x, bool new_x, Index m, Number *g)
Method to return the constraint residuals.
virtual bool eval_f(Index n, const Number *x, bool new_x, Number &obj_value)
Method to return the objective value.
Index Nt_
Number of mesh points in t direction.
Number Number Index Number Number Index Index Index Eval_F_CB eval_f
Callback function for evaluating objective function.
Number dt_
Step size in t direction.
Number lb_y_
overall lower bound on y