草庐IT

非线性最优化问题求解器Ipopt介绍

小作坊钳工 2023-04-20 原文

Ipopt(Interior Point OPTimizer)是求解大规模非线性最优化问题的求解软件。可以求解如下形式的最优化问题的(局部)最优解。
m i n ⏟ x ∈ R n     f ( x ) s . t . g L ≤ g ( x ) ≤ g U x L ≤ x ≤ x U (0) \underbrace{min}_ {x \in Rⁿ} \, \, \, f(x) \\ s.t. g_L ≤ g(x) ≤ g_U \\ x_L ≤ x ≤ x_U \tag{0} xRn minf(x)s.t.gLg(x)gUxLxxU(0)
其中, f ( x ) : R n → R f(x): R^n \rightarrow R f(x):RnR是优化目标函数, g ( x ) : R n → R m g(x):R^n \rightarrow R^m g(x):RnRm是约束函数, f ( x ) , g ( x ) f(x),g(x) f(x),g(x)可以是非线性和非凸的,但是需要是二阶微分连续的。

为了求解最优化问题,Ipopt需要更多的信息,如下:

  1. 优化问题的维度
    1. 优化变量 x x x的数目;
    2. 约束函数 g ( x ) g(x) g(x)的数目;
  2. 优化变量的边界
    1. 优化变量 x x x的边界;
    2. 约束函数 g ( x ) g(x) g(x)的边界;
  3. 优化问题的初始迭代点:
    1. 优化变量 x x x的初始值;
    2. 拉格朗日乘子的初始值(仅仅是在warm start的时候需要);
  4. 优化问题的数据结构(Structure):
    1. 约束函数 g ( x ) g(x) g(x)的雅可比矩阵的非零元素的数目;
    2. 拉格朗日函数的黑森矩阵的非零元素的数目;
    3. 约束函数 g ( x ) g(x) g(x)的雅可比稀疏矩阵的非零元素的行索引和列索引(sparsity structure,row and column indices of each of the nonzero entries );
    4. 拉格朗日函数的黑森稀疏矩阵的非零元素的行索引和列索引(sparsity structure,row and column indices of each of the nonzero entries );
  5. 优化问题函数的值:
    1. 优化目标函数 f ( x ) f(x) f(x)
    2. 优化目标函数的梯度函数 c c c;
    3. 约束函数 g ( x ) g(x) g(x)
    4. 约束函数的雅可比矩阵 ∇ g ( x ) T \nabla g(x) ^T g(x)T
    5. 拉格朗日函数的黑森矩阵 σ f ∇ 2 f ( x ) + ∑ i = 1 m λ i ∇ 2 g i ( x ) \sigma_f \nabla ^2 f(x) + \sum_{i=1}^m \lambda_i \nabla ^2 g_i(x) σf2f(x)+i=1mλi2gi(x),如果使用拟牛顿法(quasi-Newton options )则不需要此矩阵;

优化问题的维度和边界约束可以直接获得,并且来自于问题定义。初始迭代点会影响优化问题的是否收敛或者是否收敛到(局部)最优解,不同的初始值可能会导致收敛到不同的局部最优解。计算微分矩阵(雅可比矩阵和黑森矩阵)可能有一点复杂,Ipopt需要提供约束函数的雅可比矩阵和拉格朗日函数的黑森矩阵的非零元素以及他们所在的行索引和列索引,并且标准接口是下三角矩阵(黑森矩阵是对称矩阵)。矩阵的非零元素确定后,在整个求解过程中是不可变的,因此,非零元素不可以仅仅包含在初始值条件下,还需要包括在求解过程中会不为零的元素。

1. Example

f = x 1 x 4 ( x 1 + x 2 + x 3 ) + x 3 s . t .        x 1 x 2 x 3 x 4 ≥ 25 x 1 2 + x 2 2 + x 3 2 + x 4 2 = 40 1 ≤ x 1 , x 2 , x 3 , x 4 ≤ 5 x 0 = ( 1 , 5 , 5 , 1 ) (1-1) f = x_1 x_4 (x_1 + x_2 + x_3) + x_3 \\ s.t. \,\,\,\,\,\, x_1 x_2 x_3 x_4 \geq 25 \\ x^2_1 + x^2_2 + x^2_3 + x^2_4 = 40 \\ 1 \leq x_1, x_2, x_3, x_4 \leq 5 \\ x_0 = (1,5,5,1) \tag{1-1} f=x1x4(x1+x2+x3)+x3s.t.x1x2x3x425x12+x22+x32+x42=401x1,x2,x3,x45x0=(1,5,5,1)(1-1)

可以得到优化目标的梯度为:
∇ f ( x ) = [ x 1 x 4 + x 4 ( x 1 + x 2 + x 3 ) x 1 x 4 x 1 x 4 + 1 x 1 ( x 1 + x 2 + x 3 ) ] (1-2) \nabla f(x) = \begin{bmatrix} x_1 x_4 + x_4(x_1 + x_2 + x_3) \\ x_1 x_4 \\ x_1 x_4 + 1 \\ x_1 (x_1 + x_2 + x_3) \end{bmatrix} \tag{1-2} f(x)= x1x4+x4(x1+x2+x3)x1x4x1x4+1x1(x1+x2+x3) (1-2)
可以得到约束函数的雅可比矩阵为:
∇ g ( x ) = [ x 2 x 3 x 4 x 1 x 3 x 4 x 1 x 2 x 4 x 1 x 2 x 3 2 x 1 2 x 2 2 x 3 2 x 4 ] (1-3) \nabla g(x) = \begin{bmatrix} x_2 x_3 x_4 & x_1 x_3 x_4 & x_1 x_2 x_4 & x_1 x_2 x_3 \\ 2 x_1 & 2 x_2 & 2 x_3 & 2 x_4 \end{bmatrix} \tag{1-3} g(x)=[x2x3x42x1x1x3x42x2x1x2x42x3x1x2x32x4](1-3)
需要计算拉格朗日函数的黑森矩阵,如果使用拟牛顿法来近似二阶微分,则不需要提供黑森矩阵。但是黑森矩阵可以是计算有更好的鲁棒性和更快的收敛速度。NLP的拉格朗日函数定义为 f ( x ) + g ( x ) T λ f(x) + g(x)^T \lambda f(x)+g(x)Tλ,黑森矩阵定义为 ∇ 2 f ( x ) + ∑ i = 1 m λ i ∇ 2 g i ( x ) \nabla ^2 f(x) + \sum_{i=1}^m \lambda_i \nabla ^2 g_i(x) 2f(x)+i=1mλi2gi(x),然而在Ipopt中引入 σ f \sigma_f σf使Ipopt可以分别确定优化目标函数和约束函数的黑森矩阵,因此Ipopt的黑森矩阵为 σ f ∇ 2 f ( x ) + ∑ i = 1 m λ i ∇ 2 g i ( x ) \sigma_f \nabla ^2 f(x) + \sum_{i=1}^m \lambda_i \nabla ^2 g_i(x) σf2f(x)+i=1mλi2gi(x)。可以得到优化问题 ( 1 − 1 ) (1-1) (11)的黑森矩阵为:
σ f [ 2 x 4 x 4 x 4 2 x 1 + x 2 + x 3 x 4 0 0 x 1 x 4 0 0 x 1 2 x 1 + x 2 + x 3 x 1 x 1 0 ] + λ 1 [ 0 x 3 x 4 x 2 x 4 x 2 x 3 x 3 x 4 0 x 1 x 4 x 1 x 3 x 2 x 4 x 1 x 4 0 x 1 x 2 x 2 x 3 x 1 x 3 x 1 x 2 0 ] + λ 2 [ 2 0 0 0 0 2 0 0 0 0 2 0 0 0 0 2 ] (1-4) \sigma_f \begin{bmatrix} 2 x_4 & x_4 & x_4 & 2 x_1 + x_2 + x_3 \\ x_4 & 0 & 0 & x_1 \\ x_4 & 0 & 0 & x_1 \\ 2 x_1 + x_2 + x_3 & x_1 & x_1 & 0 \end{bmatrix}+ \lambda_1 \begin{bmatrix} 0 & x_3 x_4 & x_2 x_4 & x_2 x_3 \\ x_3 x_4 & 0 & x_1 x_4 & x_1 x_3 \\ x_2 x_4 & x_1 x_4 & 0 & x_1 x_2 \\ x_2 x_3 & x_1 x_3 & x_1 x_2 & 0 \end{bmatrix}+ \lambda_2 \begin{bmatrix} 2 & 0 & 0 & 0 \\ 0 & 2 & 0 & 0 \\ 0 & 0 & 2 & 0 \\ 0 & 0 & 0 & 2 \end{bmatrix} \tag{1-4} σf 2x4x4x42x1+x2+x3x400x1x400x12x1+x2+x3x1x10 +λ1 0x3x4x2x4x2x3x3x40x1x4x1x3x2x4x1x40x1x2x2x3x1x3x1x20 +λ2 2000020000200002 (1-4)
其中,第一项是优化目标函数的黑森矩阵,第二项和第三项是约束函数的黑森矩阵, λ 1 , λ 2 \lambda_1,\lambda_2 λ1,λ2是约束函数的拉格朗日乘子。

2. C++ Interface

需要继承纯虚基类Ipopt::TNLP来编写自己的求解类,并且需要重载 9 9 9Ipopt::TNLP基类的虚函数,Ipopt通过Ipopt::IpoptApplication类来求解最优化问题。

2.1 Ipopt::TNLP::get_nlp_info

   virtual bool get_nlp_info(
      Index&          n,
      Index&          m,
      Index&          nnz_jac_g,
      Index&          nnz_h_lag,
      IndexStyleEnum& index_style
   ) = 0;

Ipopt使用这个函数来确定数组的内存分配,这里如果发生问题,会引起内存泄漏等问题,很难去debug

  • n:优化变量 x x x的数目;
  • m:约束函数 g ( x ) g(x) g(x)的数目;
  • nnz_jac_g:雅可比矩阵非零元素的数目;
  • nnz_h_lag:黑森矩阵非零元素的数目;
  • index_style:稀疏矩阵的索引使用C语言风格(从 0 0 0开始),还是使用Fortran语言风格(从 1 1 1开始);

上述例子中有 4 4 4个优化变量, 2 2 2个约束函数,雅可比矩阵中的非零元素个数为 8 8 8,黑森矩阵中非零元素的个数为 16 16 16,由于是对称矩阵,因此下三角矩阵中非零元素个数为 10 10 10

// returns the size of the problem
bool HS071_NLP::get_nlp_info(
   Index&          n,
   Index&          m,
   Index&          nnz_jac_g,
   Index&          nnz_h_lag,
   IndexStyleEnum& index_style
)
{
   // The problem described in HS071_NLP.hpp has 4 variables, x[0] through x[3]
   n = 4;
 
   // one equality constraint and one inequality constraint
   m = 2;
 
   // in this example the jacobian is dense and contains 8 nonzeros
   nnz_jac_g = 8;
 
   // the Hessian is also dense and has 16 total nonzeros, but we
   // only need the lower left corner (since it is symmetric)
   nnz_h_lag = 10;
 
   // use the C style indexing (0-based)
   index_style = TNLP::C_STYLE;
 
   return true;
}

2.2 Ipopt::TNLP::get_bounds_info

   virtual bool get_bounds_info(
      Index   n,
      Number* x_l,
      Number* x_u,
      Index   m,
      Number* g_l,
      Number* g_u
   ) = 0;

Ipopt使用这个函数来确定优化变量 x x x的边界和约束函数 g ( x ) g(x) g(x)的边界。

  • n:优化变量 x x x的数目;
  • x_l:优化变量 x x x的下边界,数组;
  • x_u:优化变量 x x x的上边界,数组;
  • m:约束函数 g ( x ) g(x) g(x)的数目;
  • g_l:约束函数 g ( x ) g(x) g(x)的下边界,数组;
  • g_u:约束函数 g ( x ) g(x) g(x)的上边界,数组;

Ipopt中默认设置边界值需要在 ( − 1 0 9 , 1 0 9 ) (-10^9, 10^9) (109,109)范围内,当不在此范围时,则认为是无穷大或者无穷小。

// returns the variable bounds
bool HS071_NLP::get_bounds_info(
   Index   n,
   Number* x_l,
   Number* x_u,
   Index   m,
   Number* g_l,
   Number* g_u
)
{
   // here, the n and m we gave IPOPT in get_nlp_info are passed back to us.
   // If desired, we could assert to make sure they are what we think they are.
   assert(n == 4);
   assert(m == 2);
 
   // the variables have lower bounds of 1
   for( Index i = 0; i < 4; i++ )
   {
      x_l[i] = 1.0;
   }
 
   // the variables have upper bounds of 5
   for( Index i = 0; i < 4; i++ )
   {
      x_u[i] = 5.0;
   }
 
   // the first constraint g1 has a lower bound of 25
   g_l[0] = 25;
   // the first constraint g1 has NO upper bound, here we set it to 2e19.
   // Ipopt interprets any number greater than nlp_upper_bound_inf as
   // infinity. The default value of nlp_upper_bound_inf and nlp_lower_bound_inf
   // is 1e19 and can be changed through ipopt options.
   g_u[0] = 2e19;
 
   // the second constraint g2 is an equality constraint, so we set the
   // upper and lower bound to the same value
   g_l[1] = g_u[1] = 40.0;
 
   return true;
}	

2.3 Ipopt::TNLP::get_starting_point

   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
   ) = 0;

Ipopt使用这个函数来确定迭代优化的起点。

  • n:优化变量 x x x的数目;
  • init_x:如果是ture,则需要提供优化变量 x x x的初始值;
  • x:优化变量 x x x的初始值;

其他为dual variables的初始值,一般不用设置。在Ipopt中默认是需要设置 x x x的初始值。

// returns the initial point for the problem
bool HS071_NLP::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
)
{
   // Here, we assume we only have starting values for x, if you code
   // your own NLP, you can provide starting values for the dual variables
   // if you wish
   assert(init_x == true);
   assert(init_z == false);
   assert(init_lambda == false);
 
   // initialize to the given starting point
   x[0] = 1.0;
   x[1] = 5.0;
   x[2] = 5.0;
   x[3] = 1.0;
 
   return true;
}

2.4 Ipopt::TNLP::eval_f

   virtual bool eval_f(
      Index         n,
      const Number* x,
      bool          new_x,
      Number&       obj_value
   ) = 0;

Ipopt使用这个函数来确定优化目标函数。

  • n:优化变量 x x x的数目;
  • x:优化变量 x x x的值,用来计算 f ( x ) f(x) f(x)
  • new_x:在此之前调用的eval_*函数是否有错误发生,可以忽略;
  • obj_value f ( x ) f(x) f(x)
// returns the value of the objective function
bool HS071_NLP::eval_f(
   Index         n,
   const Number* x,
   bool          new_x,
   Number&       obj_value
)
{
   assert(n == 4);
 
   obj_value = x[0] * x[3] * (x[0] + x[1] + x[2]) + x[2];
 
   return true;
}

2.5 Ipopt::TNLP::eval_grad_f

   virtual bool eval_grad_f(
      Index         n,
      const Number* x,
      bool          new_x,
      Number*       grad_f
   ) = 0;

Ipopt使用这个函数来确定优化目标函数的梯度。

  • n:优化变量 x x x的数目;
  • x:优化变量 x x x的值,用来计算 ∇ f ( x ) \nabla f(x) f(x)
  • new_x:在此之前调用的eval_*函数是否有错误发生,可以忽略;
  • obj_value ∇ f ( x ) \nabla f(x) f(x),数组的大小和 x x x的数组大小一致;
// return the gradient of the objective function grad_{x} f(x)
bool HS071_NLP::eval_grad_f(
   Index         n,
   const Number* x,
   bool          new_x,
   Number*       grad_f
)
{
   assert(n == 4);
 
   grad_f[0] = x[0] * x[3] + x[3] * (x[0] + x[1] + x[2]);
   grad_f[1] = x[0] * x[3];
   grad_f[2] = x[0] * x[3] + 1;
   grad_f[3] = x[0] * (x[0] + x[1] + x[2]);
 
   return true;
}

2.6 Ipopt::TNLP::eval_g

   virtual bool eval_g(
      Index         n,
      const Number* x,
      bool          new_x,
      Index         m,
      Number*       g
   ) = 0;

Ipopt使用这个函数来确定约束函数 g ( x ) g(x) g(x)

  • n:优化变量 x x x的数目;
  • x:优化变量 x x x的值,用来计算 ∇ f ( x ) \nabla f(x) f(x)
  • new_x:在此之前调用的eval_*函数是否有错误发生,可以忽略;
  • m:约束函数 g ( x ) g(x) g(x)的数目;
  • g: g ( x ) g(x) g(x),数组的大小和 m m m一致;
// return the value of the constraints: g(x)
bool HS071_NLP::eval_g(
   Index         n,
   const Number* x,
   bool          new_x,
   Index         m,
   Number*       g
)
{
   assert(n == 4);
   assert(m == 2);
 
   g[0] = x[0] * x[1] * x[2] * x[3];
   g[1] = x[0] * x[0] + x[1] * x[1] + x[2] * x[2] + x[3] * x[3];
 
   return true;
}

2.7 Ipopt::TNLP::eval_jac_g

   virtual bool eval_jac_g(
      Index         n,
      const Number* x,
      bool          new_x,
      Index         m,
      Index         nele_jac,
      Index*        iRow,
      Index*        jCol,
      Number*       values
   ) = 0;

Ipopt使用这个函数来确定约束函数 g ( x ) g(x) g(x)的雅可比矩阵的非零元素的值,以及其在稀疏矩阵中的行索引值和列索引值。雅可比矩阵中的第 i i i行和第 j j j列的元素值是 g i ( x ) g_i(x) gi(x) x j x_j xj的导数。

  • n:优化变量 x x x的数目;
  • x:优化变量 x x x的值,用来计算 ∇ g ( x ) T \nabla g(x)^T g(x)T
  • new_x:在此之前调用的eval_*函数是否有错误发生,可以忽略;
  • m:约束函数 g ( x ) g(x) g(x)的数目;
  • nele_jac:雅可比矩阵非零元素的数目;
  • iRow:存储雅可比矩阵非零元素在矩阵中的行索引值,如果是C语言风格,雅可比矩阵索引值从 0 0 0开始;
  • jCol:存储雅可比矩阵非零元素在矩阵中的列索引值,如果是C语言风格,雅可比矩阵索引值从 0 0 0开始;
  • values:存储雅可比矩阵中的非零元素;

需要注意的是:①iRowjColvalues三个数组的大小是一致的,并且其储存的值应该和雅可比矩阵非零元素的行索引值、列索引值和非零元素值相对应;②数组iRowjCol只需要被填写一次,即第一次调用此函数时填写iRowjCol,第一次调用时xvalues都是null,当Ipopt需要values的值时,传递iRowjCol将会是null,此时对values的值进行填写。

// return the structure or values of the Jacobian
bool HS071_NLP::eval_jac_g(
   Index         n,
   const Number* x,
   bool          new_x,
   Index         m,
   Index         nele_jac,
   Index*        iRow,
   Index*        jCol,
   Number*       values
)
{
   assert(n == 4);
   assert(m == 2);
 
   if( values == NULL )
   {
      // return the structure of the Jacobian
 
      // this particular Jacobian is dense
      iRow[0] = 0;
      jCol[0] = 0;
      iRow[1] = 0;
      jCol[1] = 1;
      iRow[2] = 0;
      jCol[2] = 2;
      iRow[3] = 0;
      jCol[3] = 3;
      iRow[4] = 1;
      jCol[4] = 0;
      iRow[5] = 1;
      jCol[5] = 1;
      iRow[6] = 1;
      jCol[6] = 2;
      iRow[7] = 1;
      jCol[7] = 3;
   }
   else
   {
      // return the values of the Jacobian of the constraints
 
      values[0] = x[1] * x[2] * x[3]; // 0,0
      values[1] = x[0] * x[2] * x[3]; // 0,1
      values[2] = x[0] * x[1] * x[3]; // 0,2
      values[3] = x[0] * x[1] * x[2]; // 0,3
 
      values[4] = 2 * x[0]; // 1,0
      values[5] = 2 * x[1]; // 1,1
      values[6] = 2 * x[2]; // 1,2
      values[7] = 2 * x[3]; // 1,3
   }
 
   return true;
}

2.8 Ipopt::TNLP::eval_h

   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
   )

Ipopt使用这个函数来确定拉格朗日函数黑森矩阵的非零元素的值,以及其在稀疏矩阵中的行索引值和列索引值。

  • n:优化变量 x x x的数目;
  • x:优化变量 x x x的值,用来计算 ∇ g ( x ) T \nabla g(x)^T g(x)T
  • new_x:在此之前调用的eval_*函数是否有错误发生,可以忽略;
  • obj_factor σ f \sigma_f σf
  • m:约束函数 g ( x ) g(x) g(x)的数目;
  • lambda:拉格朗日乘子 λ \lambda λ
  • new_lambda:如果之前调用的函数使用相同的 λ \lambda λ则为false,一般忽略;
  • nele_hess:黑森矩阵非零元素的个数(下三角矩阵);
  • iRow:存储黑森矩阵非零元素在矩阵中的行索引值,如果是C语言风格,黑森矩阵索引值从 0 0 0开始;
  • jCol:存储黑森矩阵非零元素在矩阵中的列索引值,如果是C语言风格,黑森矩阵索引值从 0 0 0开始;
  • values:存储黑森矩阵中的非零元素的值;

需要注意的是:①iRowjColvalues三个数组的大小是一致的,并且其储存的值应该和黑森矩阵非零元素的行索引值、列索引值和非零元素值相对应;②数组iRowjCol只需要被填写一次,即第一次调用此函数时填写iRowjCol,第一次调用时xlambdavalues都是null,当Ipopt需要values的值时,传递iRowjCol将会是null,此时对values的值进行填写;③由于黑森矩阵是对称阵,Ipopt使用下三角矩阵;④Ipopt默认是需要黑森矩阵的,当使用拟牛顿法时,则不需要黑森矩阵。

在此例中,黑森矩阵是稠密的,但是仍然使用稀疏矩阵来表示。

//return the structure or values of the Hessian
bool HS071_NLP::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
)
{
   assert(n == 4);
   assert(m == 2);
 
   if( values == NULL )
   {
      // return the structure. This is a symmetric matrix, fill the lower left
      // triangle only.
 
      // the hessian for this problem is actually dense
      Index idx = 0;
      for( Index row = 0; row < 4; row++ )
      {
         for( Index col = 0; col <= row; col++ )
         {
            iRow[idx] = row;
            jCol[idx] = col;
            idx++;
         }
      }
 
      assert(idx == nele_hess);
   }
   else
   {
      // return the values. This is a symmetric matrix, fill the lower left
      // triangle only
 
      // fill the objective portion
      values[0] = obj_factor * (2 * x[3]); // 0,0
 
      values[1] = obj_factor * (x[3]);     // 1,0
      values[2] = 0.;                      // 1,1
 
      values[3] = obj_factor * (x[3]);     // 2,0
      values[4] = 0.;                      // 2,1
      values[5] = 0.;                      // 2,2
 
      values[6] = obj_factor * (2 * x[0] + x[1] + x[2]); // 3,0
      values[7] = obj_factor * (x[0]);                   // 3,1
      values[8] = obj_factor * (x[0]);                   // 3,2
      values[9] = 0.;                                    // 3,3
 
      // add the portion for the first constraint
      values[1] += lambda[0] * (x[2] * x[3]); // 1,0
 
      values[3] += lambda[0] * (x[1] * x[3]); // 2,0
      values[4] += lambda[0] * (x[0] * x[3]); // 2,1
 
      values[6] += lambda[0] * (x[1] * x[2]); // 3,0
      values[7] += lambda[0] * (x[0] * x[2]); // 3,1
      values[8] += lambda[0] * (x[0] * x[1]); // 3,2
 
      // add the portion for the second constraint
      values[0] += lambda[1] * 2; // 0,0
 
      values[2] += lambda[1] * 2; // 1,1
 
      values[5] += lambda[1] * 2; // 2,2
 
      values[9] += lambda[1] * 2; // 3,3
   }
 
   return true;
}

2.9 Ipopt::TNLP::finalize_solution

   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
   ) = 0;

Ipopt使用这个函数来得到最优化问题的求解结果,对其重要的值进行介绍。

  • status:求解器的状态;
    • SUCCESS:在满足收敛条件的情况下,找到局部最优解;
    • MAXITER_EXCEEDED:超出最大迭代次数;
    • CPUTIME_EXCEEDED:超出最大求解时间;
    • STOP_AT_ACCEPTABLE_POINT:求解收敛在某点,不满足期望的容差,但是在可接受范围内;
    • LOCAL_INFEASIBILITY:在可行域内找不到最优解,一般是由于bounds和约束设置不合理导致的;
  • x:优化变量 x x x的局部最优解的值;
void HS071_NLP::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
)
{
   // here is where we would store the solution to variables, or write to a file, etc
   // so we could use the solution.
 
   // For this example, we write the solution to the console
   std::cout << std::endl << std::endl << "Solution of the primal variables, x" << std::endl;
   for( Index i = 0; i < n; i++ )
   {
      std::cout << "x[" << i << "] = " << x[i] << std::endl;
   }
 
   std::cout << std::endl << std::endl << "Solution of the bound multipliers, z_L and z_U" << std::endl;
   for( Index i = 0; i < n; i++ )
   {
      std::cout << "z_L[" << i << "] = " << z_L[i] << std::endl;
   }
   for( Index i = 0; i < n; i++ )
   {
      std::cout << "z_U[" << i << "] = " << z_U[i] << std::endl;
   }
 
   std::cout << std::endl << std::endl << "Objective value" << std::endl;
   std::cout << "f(x*) = " << obj_value << std::endl;
 
   std::cout << std::endl << "Final value of the constraints:" << std::endl;
   for( Index i = 0; i < m; i++ )
   {
      std::cout << "g(" << i << ") = " << g[i] << std::endl;
   }
}

2.10 main function

上述对Ipopt::TNLP的函数进行了重载,但是需要编写调用Ipopt的函数来执行求解。

#include "IpIpoptApplication.hpp"
#include "hs071_nlp.hpp"
 
#include <iostream>
 
using namespace Ipopt;
 
int main(
   int    /*argv*/,
   char** /*argc*/
)
{
   // Create a new instance of your nlp
   //  (use a SmartPtr, not raw)
   SmartPtr<TNLP> mynlp = new HS071_NLP();
 
   // Create a new instance of IpoptApplication
   //  (use a SmartPtr, not raw)
   // We are using the factory, since this allows us to compile this
   // example with an Ipopt Windows DLL
   SmartPtr<IpoptApplication> app = IpoptApplicationFactory();
 
   // Change some options
   // Note: The following choices are only examples, they might not be
   //       suitable for your optimization problem.
   app->Options()->SetNumericValue("tol", 3.82e-6);
   app->Options()->SetStringValue("mu_strategy", "adaptive");
   app->Options()->SetStringValue("output_file", "ipopt.out");
   // The following overwrites the default name (ipopt.opt) of the options file
   // app->Options()->SetStringValue("option_file_name", "hs071.opt");
 
   // Initialize the IpoptApplication and process the options
   ApplicationReturnStatus status;
   status = app->Initialize();
   if( status != Solve_Succeeded )
   {
      std::cout << std::endl << std::endl << "*** Error during initialization!" << std::endl;
      return (int) status;
   }
 
   // Ask Ipopt to solve the problem
   status = app->OptimizeTNLP(mynlp);
 
   if( status == Solve_Succeeded )
   {
      std::cout << std::endl << std::endl << "*** The problem solved!" << std::endl;
   }
   else
   {
      std::cout << std::endl << std::endl << "*** The problem FAILED!" << std::endl;
   }
 
   // As the SmartPtrs go out of scope, the reference count
   // will be decremented and the objects will automatically
   // be deleted.
 
   return (int) status;
}

有关非线性最优化问题求解器Ipopt介绍的更多相关文章

  1. ruby - 在 64 位 Snow Leopard 上使用 rvm、postgres 9.0、ruby 1.9.2-p136 安装 pg gem 时出现问题 - 2

    我想为Heroku构建一个Rails3应用程序。他们使用Postgres作为他们的数据库,所以我通过MacPorts安装了postgres9.0。现在我需要一个postgresgem并且共识是出于性能原因你想要pggem。但是我对我得到的错误感到非常困惑当我尝试在rvm下通过geminstall安装pg时。我已经非常明确地指定了所有postgres目录的位置可以找到但仍然无法完成安装:$envARCHFLAGS='-archx86_64'geminstallpg--\--with-pg-config=/opt/local/var/db/postgresql90/defaultdb/po

  2. ruby - 通过 rvm 升级 ruby​​gems 的问题 - 2

    尝试通过RVM将RubyGems升级到版本1.8.10并出现此错误:$rvmrubygemslatestRemovingoldRubygemsfiles...Installingrubygems-1.8.10forruby-1.9.2-p180...ERROR:Errorrunning'GEM_PATH="/Users/foo/.rvm/gems/ruby-1.9.2-p180:/Users/foo/.rvm/gems/ruby-1.9.2-p180@global:/Users/foo/.rvm/gems/ruby-1.9.2-p180:/Users/foo/.rvm/gems/rub

  3. ruby - 通过 RVM (OSX Mountain Lion) 安装 Ruby 2.0.0-p247 时遇到问题 - 2

    我的最终目标是安装当前版本的RubyonRails。我在OSXMountainLion上运行。到目前为止,这是我的过程:已安装的RVM$\curl-Lhttps://get.rvm.io|bash-sstable检查已知(我假设已批准)安装$rvmlistknown我看到当前的稳定版本可用[ruby-]2.0.0[-p247]输入命令安装$rvminstall2.0.0-p247注意:我也试过这些安装命令$rvminstallruby-2.0.0-p247$rvminstallruby=2.0.0-p247我很快就无处可去了。结果:$rvminstall2.0.0-p247Search

  4. ruby - Fast-stemmer 安装问题 - 2

    由于fast-stemmer的问题,我很难安装我想要的任何ruby​​gem。我把我得到的错误放在下面。Buildingnativeextensions.Thiscouldtakeawhile...ERROR:Errorinstallingfast-stemmer:ERROR:Failedtobuildgemnativeextension./System/Library/Frameworks/Ruby.framework/Versions/2.0/usr/bin/rubyextconf.rbcreatingMakefilemake"DESTDIR="cleanmake"DESTDIR=

  5. ruby - 安装 Ruby 时遇到问题(无法下载资源 "readline--patch") - 2

    当我尝试安装Ruby时遇到此错误。我试过查看this和this但无济于事➜~brewinstallrubyWarning:YouareusingOSX10.12.Wedonotprovidesupportforthispre-releaseversion.Youmayencounterbuildfailuresorotherbreakages.Pleasecreatepull-requestsinsteadoffilingissues.==>Installingdependenciesforruby:readline,libyaml,makedepend==>Installingrub

  6. java - 从 JRuby 调用 Java 类的问题 - 2

    我正在尝试使用boilerpipe来自JRuby。我看过guide从JRuby调用Java,并成功地将它与另一个Java包一起使用,但无法弄清楚为什么同样的东西不能用于boilerpipe。我正在尝试基本上从JRuby中执行与此Java等效的操作:URLurl=newURL("http://www.example.com/some-location/index.html");Stringtext=ArticleExtractor.INSTANCE.getText(url);在JRuby中试过这个:require'java'url=java.net.URL.new("http://www

  7. ruby-on-rails - 简单的 Ruby on Rails 问题——如何将评论附加到用户和文章? - 2

    我意识到这可能是一个非常基本的问题,但我现在已经花了几天时间回过头来解决这个问题,但出于某种原因,Google就是没有帮助我。(我认为部分问题在于我是一个初学者,我不知道该问什么......)我也看过O'Reilly的RubyCookbook和RailsAPI,但我仍然停留在这个问题上.我找到了一些关于多态关系的信息,但它似乎不是我需要的(尽管如果我错了请告诉我)。我正在尝试调整MichaelHartl'stutorial创建一个包含用户、文章和评论的博客应用程序(不使用脚手架)。我希望评论既属于用户又属于文章。我的主要问题是:我不知道如何将当前文章的ID放入评论Controller。

  8. 【高数】用拉格朗日中值定理解决极限问题 - 2

    首先回顾一下拉格朗日定理的内容:函数f(x)是在闭区间[a,b]上连续、开区间(a,b)上可导的函数,那么至少存在一个,使得:通过这个表达式我们可以知道,f(x)是函数的主体,a和b可以看作是主体函数f(x)中所取的两个值。那么可以有,  也就意味着我们可以用来替换 这种替换可以用在求某些多项式差的极限中。方法: 外层函数f(x)是一致的,并且h(x)和g(x)是等价无穷小。此时,利用拉格朗日定理,将原式替换为 ,再进行求解,往往会省去复合函数求极限的很多麻烦。使用要注意:1.要先找到主体函数f(x),即外层函数必须相同。2.f(x)找到后,复合部分是等价无穷小。3.要满足作差的形式。如果是加

  9. Unity 热更新技术 | (三) Lua语言基本介绍及下载安装 - 2

    ?博客主页:https://xiaoy.blog.csdn.net?本文由呆呆敲代码的小Y原创,首发于CSDN??学习专栏推荐:Unity系统学习专栏?游戏制作专栏推荐:游戏制作?Unity实战100例专栏推荐:Unity实战100例教程?欢迎点赞?收藏⭐留言?如有错误敬请指正!?未来很长,值得我们全力奔赴更美好的生活✨------------------❤️分割线❤️-------------------------

  10. SPI接收数据异常问题总结 - 2

    SPI接收数据左移一位问题目录SPI接收数据左移一位问题一、问题描述二、问题分析三、探究原理四、经验总结最近在工作在学习调试SPI的过程中遇到一个问题——接收数据整体向左移了一位(1bit)。SPI数据收发是数据交换,因此接收数据时从第二个字节开始才是有效数据,也就是数据整体向右移一个字节(1byte)。请教前辈之后也没有得到解决,通过在网上查阅前人经验终于解决问题,所以写一个避坑经验总结。实际背景:MCU与一款芯片使用spi通信,MCU作为主机,芯片作为从机。这款芯片采用的是它规定的六线SPI,多了两根线:RDY和INT,这样从机就可以主动请求主机给主机发送数据了。一、问题描述根据从机芯片手

随机推荐