非线性优化--NLopt算法使用及C++实例

看这篇之前建议先看这篇，里面讲了非线性优化的原理即相关名词的概念，然后介绍了NLopt的使用方法，这个方法是基于C语言的，本片介绍一个NLopt的实例，用的C++语言。

在实例之前，先介绍下NLopt支持的算法，以及算法使用的注意事项

NLopt 支持的算法

NLopt 包含很多种不同的优化算法。

在头文件里面算法名称的枚举类型为

  enum algorithm {
     GN_DIRECT = 0,
     GN_DIRECT_L,
     GN_DIRECT_L_RAND,
     GN_DIRECT_NOSCAL,
     GN_DIRECT_L_NOSCAL,
     GN_DIRECT_L_RAND_NOSCAL,
     GN_ORIG_DIRECT,
     GN_ORIG_DIRECT_L,
     GD_STOGO,
     GD_STOGO_RAND,
     LD_LBFGS_NOCEDAL,
     LD_LBFGS,
     LN_PRAXIS,
     LD_VAR1,
     LD_VAR2,
     LD_TNEWTON,
     LD_TNEWTON_RESTART,
     LD_TNEWTON_PRECOND,
     LD_TNEWTON_PRECOND_RESTART,
     GN_CRS2_LM,
     GN_MLSL,
     GD_MLSL,
     GN_MLSL_LDS,
     GD_MLSL_LDS,
     LD_MMA,
     LN_COBYLA,
     LN_NEWUOA,
     LN_NEWUOA_BOUND,
     LN_NELDERMEAD,
     LN_SBPLX,
     LN_AUGLAG,
     LD_AUGLAG,
     LN_AUGLAG_EQ,
     LD_AUGLAG_EQ,
     LN_BOBYQA,
     GN_ISRES,
     AUGLAG,
     AUGLAG_EQ,
     G_MLSL,
     G_MLSL_LDS,
     LD_SLSQP,
     LD_CCSAQ,
     GN_ESCH,
     NUM_ALGORITHMS /*不是一种算法 只是算法的数量*/
  };

命名规律：

G/L代表的就是 全局（global）或者局部（local）优化
N/D代表的就是 不需导数 或者 需要导数 的优化

例如 LN_COBYLA 就是用的 COBYLA 算法 ，然后该算法用于局部（L）无需导数（N）的优化

算法选择

那么我们在使用的时候用哪种算法呢？随意选择吗？NO
对于一个确定的数学模型的优化问题，比较好的方式是 通过比较几种可用的算法，再确定选哪种，因为没有最佳算法，不同优化问题对于最佳的解决方式不同。

在比较算法的时候也要注意，不同算法适配的 函数值容差和参数容差 也应该是不同的。

有一种比较好的比较两种算对于一个数学模型的优劣的方式：
先运行一个算法得到最小值，然后运行第二个算法的时候设定就收敛到那个值，然后比较运行的时间

选择全局优化要注意的问题

当前，所有全局优化算法都要求对所有优化参数指定边界约束。
并且一定要注意的是
在这些算法中 只有 ISRES, AGS, 和 ORIG_DIRECT 支持非线性不等式约束
只有 ISRES 支持非线性等式约束

有一种很好的做法就是通过全局算法找到一个最优点，然后以这个最优点作为局部优化的起点，使结果更准确

全局优化算法在找到最优点时 比局部优化算法要 花的精力多很多

Code

数学模型：

求 ：max（lnx1+lnx2）    目标函数
约束： p1x1+p2x2=5        等式约束
            x1<=x2                  不等式约束                
            x1>=0;x2>=0         边界约束

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

#include <nlopt.hpp> //nlopt的头文件

引用头文件
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

//声明一个 是否使用梯度的 全局变量
bool grad_bool=1;

声明一个 是否使用梯度的 全局变量
这样方便在切换 基于梯度的算法和不需要导数的算法的时候可以 改下这个变量就可以
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

double myfunc(const std::vector<double>& x,std::vector<double>&grad,void* f_data )
{
    double result;//声明结果
    if (grad_bool)
    {
        grad[0] = 1 / x[0];
        grad[1] = 1 / x[1];
    }
    result = log(x[0])+log(x[1]);//计算结果
    return result;//返回结果
}

目标函数是 max（lnx1+x2）
定义目标函数 注意函数的参数的形式 是不可以变化的
参数x 是 要优化的参数向量
grad返回为此时最优化参数的梯度 就是数学模型里里 对x求偏导的结果 如果要求偏导那么 lnx1+lnx2 对x1求偏导，及1/x1 也就是 代码的 grad[0] = 1 / x[0]; 这个地方
f_data 是要传入的参数
返回值 为在一个x向量下的目标函数的值
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

double myconstraint(const std::vector<double>& x,std::vector<double>&grad,void* f_data)
{
    double *p = (double *) f_data;
    if (grad_bool)
    {
        grad[0] = p[0];
        grad[1] = p[1];
    }
    double result;//声明结果    
    result = x[0] * p[0] + x[1] * p[1] - 5;//计算结果
    return result;//返回结果
}

定义等式约束函数
数学模型的等式是：p1x1+p2x2=5
p1和p2以外部参数的方式传进来

最好的方式换种写法，这样什么样的参数格式都可以套用

typedef struct{
    double p1,p2;
}my_constraint_data;

声明要传入参数的结构体

double myconstraint(const std::vector<double>& x,std::vector<double>&grad,void* f_data)
{
    my_constraint_data *d = (my_constraint_data*)f_data;
    //double *p = (double *) f_data;
    double p1 = d->p1;
    double p2 = d->p2;
    if (grad_bool)
    {
        grad[0] = p1;
        grad[1] = p2;
    }
    double result;//声明结果    
    result = x[0] *  p1  + x[1] * p2 - 5;//计算结果
    return result;//返回结果
}

然后取参数的时候这样就可以

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

double myinconstraint(const  std::vector<double>& x,std::vector<double>&grad,void* f_data)
{
    if (grad_bool)
    {
        grad[0] = 1;
        grad[1] = -1;
    }

    double result;//声明结果    
    result = x[0]  - x[1] ;//计算结果
    return result;//返回结果
}

不等式的约束函数
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

int main(int argc, char** argv) 
{
    //初始化节点
    ros::init(argc, argv, "lidar_align");

    //声明两个句柄   
    ros::NodeHandle nh;

然后就是main函数的，我是在ros下做的，所以加个节点的初始化
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    //数学模型里的参数
    //double p[2]={1,2};
    my_constraint_data data = {1,2};

    //待优化 求取的参数
    std::vector<double> x{1,1};

声明数学模型里面 等式约束的 参数 my_constraint_data这个结构体在上面定义了
声明 待优化 求取的参数 1,1相当于初始化了
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    //声明一个优化器
    nlopt::opt opter;

声明一个NLopt的优化器

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /*局部优化算法*/
    //opter = nlopt::opt(nlopt::LD_SLSQP,2);  //有导数
    //opter = nlopt::opt(nlopt::LN_COBYLA,2); //无导数
    /*全局优化算法*/
    opter = nlopt::opt(nlopt::GN_ISRES,2);    //无导数  是可以的   x1 = 1.66667    x2= 1.66667  fmax = 1.02165

给优化器设置使用什么优化算法，并设置因子个数，就是x的个数
优化算法在nlopt.hpp里面有， 可以从里面选择想用的

但是注意并不是所有算法都可以计算出来，例如这个数学模型有 不等式约束和等式约束的 全局优化算法里面只有ISRES可以使用

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /*优化参数的边界约束*/
    std::vector<double> lb {1,1};//注意参数的个数要对应上
    std::vector<double> ub {10000,10000};//注意参数的个数要对应上
        //设置 参数 边界
    opter.set_lower_bounds(lb);//设置参数下限
    opter.set_upper_bounds(ub);//设置参数上限

优化参数的边界约束

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /*设置优化精度*/
    double tol = 1e-8;//容差
    opter.set_xtol_rel(tol);
    opter.set_force_stop(tol);

设置优化精度
控制优化什么时候停止
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /*设置目标函数  其中目标函数是自己定义的 */
    // myfunc 是自己定义的函数名字 ，第二个参数为 可传入数据 没有数据则为NULL
    opter.set_max_objective(myfunc,NULL);

设置目标函数

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /*设置数学模型的 等式 约束 和 不等式 约束   */
    //myconstraint 是自己定义的函数名字  ,第二个参数为 可传入数据,
    opter.add_equality_constraint(myconstraint,&data,tol);
    opter.add_inequality_constraint(myinconstraint, NULL,tol);

设置数学模型的 等式 约束 和 不等式 约束
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

    /* 获取结果 最大值储存在f_max内 对应的向量储存在x内*/
    double f_max = -10000;
    nlopt::result res = opter.optimize(x,f_max);

执行优化计算
～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

 std::cout<<"x1 = " <<x[0]<<"    x2= "<<x[1]<<"  fmax = "<<f_max<<std::endl;

打印计算结果

～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～～

Result