高斯伪谱C++封装库开源！

Windows x64/86 C++无依赖运行高斯伪谱法求解最优控制问题，你只需要ElegantGP!

Author: Y. F. Zhang His Github: https://github.com/ZYunfeii

写在前面

这个库在你下载它的那一时刻起不再依赖任何其他代码，直接可用来构建C++的最优控制问题并进行求解。我还写了一个visual studio使用该库的demo项目，供学习。
项目主要基于Lpopc进行封装，编译不易，下载地址：https://download.csdn.net/download/weixin_43145941/88817667

文件简述

arma: 矩阵第三方库 https://gitlab.com/conradsnicta/armadillo-code

Lpopc: 高斯伪谱法库 https://sourceforge.net/projects/lpopc/

Debug_win64: Debug版本64位库

Release_win64: Release版本64位库

MKL: MKL相关库和一个intel导入库(libiomp5md.lib)

对于库中文件解释：

Ipopt-vc8.dll：ipopt动态库

Ipopt-vc8.lib：ipopt导入库

liblpopc.lib：高斯伪谱封装库

使用

库的介绍

../Example文件夹中给出了一个经典的轨迹优化案例visual studio项目。

使用C++版本解决高斯伪谱问题需要的库有（Release）：

Ipopt-vc8.lib
liblpopc.lib
mkl_intel_lp64.lib
mkl_intel_thread.lib
mkl_core.lib
libiomp5md.lib

第1个库为ipopt库（由https://github.com/coin-or/Ipopt 编译）。第2个库为高斯伪谱库（由 https://sourceforge.net/projects/lpopc/ 编译）。第3，4，5个库是MKL的静态库，这里我直接将其拷贝过来了，无需使用者自行安装。缺点就是这几个静态库十分臃肿。第6个库是intel相关库，我也直接拷贝过来了。

MKL全称Intel Math Kernel Library，是由Intel 公司开发的，专门用于矩阵计算的库。

visual studio项目配置

遵循Debug对应Base下Debug库，Release对应Base下Release库，编译平台选x64。
VC++目录>>包含目录：

$(SolutionDir)..\ElegantGP\Lpopc\Common
$(SolutionDir)..\ElegantGP\Lpopc\Core
$(SolutionDir)..\ElegantGP\Lpopc\SparseMatrix
$(SolutionDir)..\ElegantGP\arma\include

具体路径根据用户Base位置确定。

$(SolutionDir)..\ElegantGP\Debug_win64
$(SolutionDir)..\ElegantGP\MKL

具体路径根据用户库位置确定。再次强调，Debug和Release需对应。

链接器>>输入>>附加依赖项：

Ipopt-vc8.lib
liblpopc.lib
mkl_intel_lp64.lib
mkl_intel_thread.lib
mkl_core.lib
libiomp5md.lib

输入上述库名称。

C/C++>>代码生成>>运行库：选择多线程调试(/MTd)

其余可根据用户需求进一步细优化配置，至此，可进行项目生成。

将Ipopt-vc8.dll和libiomp5md.dll拷贝到生成的可执行文件夹下（也可以把dll配置到环境变量，但保险起见使用前者可以保证本库的dll最先被找到，因为可能你的电脑上存在同名dll，据我所知，matlab里面有，如果你配了matlab的环境变量，那很可能找到它的同名dll导致exe运行时出现程序定位点错误）【重要】
双击执行生成的example.exe

[可选] C/C++>>优化选用最大优化（优化速度）

库函数

高斯伪谱mesh refine方法选用hp-Liu（hp方法求解效率不行）:

app->Options()->SetStringValue("mesh-refine-methods", "hp-Liu");

最大网格数设置：

app->Options()->SetIntegerValue("max-grid-num", 120);

误差设置：

app->Options()->SetNumericValue("finite-difference-tol", 1e-3);
app->Options()->SetNumericValue("desired-relative-error", 1e-3);

求解结果

txt形式

在exe文件目录下生成state time control文件，其为轨迹优化结果。

代码形式

我对原库进行了修改，可直接从应用层获取求解结果：

app->algorithm_->cd_data_->result[0].get()->state; // mat形式的state，0表示phase编号，从0开始

具体可获取的结果见如下结构体：

struct SolutionData
{
    vec time;
    mat state;
    mat control;
    mat parameter;
    mat costate;
    mat pathmult;
    mat Hamiltonian;
    double mayerCost;
    double lagrangeCost;

};

lagrangeCost表示积分型代价函数值，mayerCost就是传统意义上不带积分的代价函数值。

性能相关

cpu: Intel i7-11700 16核

HyperSensitive轨迹优化Release版本求解0.581s（Release进行了编译运算优化）。
HyperSensitive轨迹优化Debug版本求解0.967s。
HyperSensitive轨迹优化Matlab2019b相同初始值求解7.937s。

请添加图片描述

计算误差在1e-10级别。

另一个例子：

高超声速飞行器再入轨迹优化问题：Matlab2019b求解35s，ElegantGP只需要1.7s。

不足

liblpopc.lib比较臃肿导致编译出来的可执行文件达几十MBytes。
MKL的库也十分臃肿，为了打包我都拷贝过来了，但是文件过大。

关于从头编译lpopc库

Linux下lpopc库的编译还是较为容易的，但也不是非常的容易。而Windows下该项目的编译可以用困难重重形容。

Lpopc的作者在文档LpopcDoc.pdf（\lpopc-master\Lpopc\doc）中给出了其编译流程，但仍旧有许多不一致。感兴趣的读者可以自行尝试编译：

git下来Ipopt的项目，进入Ipopt-3.12.3\Ipopt\MSVisualStudio\v8-ifort，最终是要把Ipopt-vc8项目编译出来。但它的编译依赖解决方案中CoinMetis，CoinMumpsC，CoinMumpsF90，IpOptFor项目编译出来的静态库，因此需要先编译这几个项目。
1中提到的需要先编译的项目中有Fortan项目，这里需要安装Fortran编译器ifort，除此之外MKL库也是必须的。
CoinMetis，CoinMumpsC，CoinMumpsF90的编译依赖METIS和MUMPS项目，需要下载源码，把1中项目.F文件用源码替代。
使用ipopt编译好的库对liblpopc库进行编译，这里liblpopc的vs项目源文件少添加了关于hpLiu的mesh方法cpp和hpp文件。
使用liblpopc的库对高斯伪谱优化问题编译。

这里只是非常简略得叙述了下编译过程，实际上有很多细微的问题，不再赘述。

对原库的细节修改

Lpopc原本是求解完最优问题后通过arma的接口将结果写入磁盘，这不利于将GP嵌入自己的算法作为中间环节。因此，我将LpLpopcAlgorithm.hpp文件中LpopcAlgorithm类的私有变量cd_data_改为public，同时将LpLpopcApplication.hpp文件中algorithm_改为public。这一改动不够优雅但无伤大雅。