使用Armadillo和qpOASES时抛出异常

Exception throw when using Armadillo and qpOASES

本文关键字：抛出异常 qpOASES Armadillo 使用      更新时间：2023-10-16

我有一个二次规划优化问题，我正在用qpOASES解决。在那里存在一个矩阵X，我需要先决条件，所以我使用Armadillo和常规arma::pinv从那里为了计算moore - penrose伪逆。

问题:我把矩阵X写在一个文件中，然后在一个不依赖于qpOASES的单独程序(比如test.cpp)中读取它。常规pinv运行良好。

#include <iostream> 
#include <fstream>
#include <armadillo>                                                  
#include <string>
using namespace std;                                       
using namespace arma;                                      
int main(){

    // Read design matrix.
    int NRows = 199;
    int NFields = 26;
    string flname_in = "chol_out_2_data"; 
    mat  A (NRows,NFields);                              
    for (int i=0; i < NRows; ++i) 
         for (int j=0; j < NFields; ++j)
                    myin >> A(i,j) ;
    // Calculate pseudoinverse  
    mat M;
    pinv(M,A);    // <========= THIS fails when I use flag: -lqpOASES
 }

当我在执行QP优化的文件中包含相同的例程时(例如true_QP.cpp)，我得到一个运行时错误，因为pinv无法计算伪逆。我已经做了大量的测试，文件在OK中读取，值是相同的。

我以以下方式跟踪了冲突的问题:我编译了不依赖于qpOASES (test.cpp -如上所述)的程序，也使用标志- lqpoases，然后，代码给出运行时错误。

也就是编译:

g++ test.cpp -o test.xxx -larmadillo

运行好:

./test.xxx
编译:

g++ test.cpp -o test.xxx -larmadillo -lqpOASES

抛出异常(由于计算pinv失败):

./test.xxx

因此我怀疑有一些冲突-似乎使用- lqpoases也会影响犰狳的一些标志?什么好主意吗?是否有一些LAPACK/BLAS的依赖或一些内部标志可能会改变犰狳的设置?谢谢你的宝贵时间。

下面是arma::pinv函数的文档:http://arma.sourceforge.net/docs.html pinv

我已经通过从Eigen而不是Armadillo计算pinv解决了这个问题。

基于这个bug报告，我为Eigen使用的函数定义:http://eigen.tuxfamily.org/bz/show_bug.cgi?id=257

: 

template<typename _Matrix_Type_>
Eigen::MatrixXd pinv(const _Matrix_Type_ &a, double epsilon =std::numeric_limits<double>::epsilon())
    {
    Eigen::JacobiSVD< _Matrix_Type_ > svd(a ,Eigen::ComputeThinU | Eigen::ComputeThinV);
    double tolerance = epsilon * std::max(a.cols(), a.rows()) *svd.singularValues().array().abs()(0);
    return
            svd.matrixV() *  (svd.singularValues().array().abs() > tolerance).select(svd.singularValues().array().inverse(), 0).matrix().asDiagonal() * svd.matrixU().adjoint();
}
相关文章：