语言, C/C++专题

C++程序的OpenBLAS链接和编译

这是之前相关的几篇:

本篇介绍 C++ 代码实现 OpenBLAS 的链接和编译。

一、C++ 直接调用 OpenBLAS

C++ 代码调用 OpenBLAS 矩阵求逆的例子 openblas_test.cpp:

#include <iostream>
#include <vector>
#include <cmath>

extern "C" {
    void dgetrf_(int* m, int* n, double* a, int* lda, int* ipiv, int* info);
    void dgetri_(int* n, double* a, int* lda, int* ipiv, double* work, int* lwork, int* info);
}

int main() {
    const int N = 3;
    
    double A_data[] = {
        1, 0, 5,   // 第1列
        2, 1, 6,   // 第2列
        3, 4, 0    // 第3列
    };
    std::vector<double> A(A_data, A_data + 9);

    std::cout << "Original matrix (column-major):\n";
    for (int i = 0; i < N; ++i) {
        for (int j = 0; j < N; ++j) {
            std::cout << A[j*N + i] << " ";
        }
        std::cout << "\n";
    }

    // LU 分解
    std::vector<int> ipiv(N);
    int info;
    int n = N;
    int lda = N;
    
    dgetrf_(&n, &n, &A[0], &lda, &ipiv[0], &info);
    if (info != 0) {
        std::cerr << "LU decomposition failed! (info=" << info << ")\n";
        return 1;
    }

    // 计算逆矩阵
    std::vector<double> work(N);
    int lwork = N; // 关键:必须是可变变量
    dgetri_(&n, &A[0], &lda, &ipiv[0], &work[0], &lwork, &info);
    if (info != 0) {
        std::cerr << "Matrix inversion failed! (info=" << info << ")\n";
        return 1;
    }

    std::cout << "\nInverse matrix:\n";
    for (int i = 0; i < N; ++i) {
        for (int j = 0; j < N; ++j) {
            std::cout << A[j*N + i] << " ";
        }
        std::cout << "\n";
    }

    return 0;
}

编译命令和 Fortran 的命令差不多,参考:Fortran程序的OpenBLAS链接和编译

默认 OpenBLAS 环境,动态链接编译:

g++ -o test openblas_test.cpp -lopenblas

指定 OpenBLAS 路径,静态链接 OpenBLAS,动态链接 pthread 库,编译:

g++ -o test openblas_test.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread

指定 OpenBLAS 路径,完全静态链接编译:

g++ -static -o test openblas_test.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread

可以查看编译后的可执行文件 test 的大小,静态链接编译的文件会更大些:

32K     test_1
316K    test_2
2.8M    test_3

编译后运行:

./test

运行结果:

Original matrix (column-major):
1 2 3
0 1 4
5 6 0

Inverse matrix:
-24 18 5
20 -15 -4
-5 4 1

二、通过 Eigen 调用 OpenBLAS(推荐)

这里通过简单修改代码,通过 Eigen 实现调用 OpenBLAS,而不是使用 Eigen 自带的默认 BLAS。Eigen 环境的安装以及原代码示例参考这篇:C++线性代数库Eigen的下载和代码示例

代码例子 openblas_test_with_eigen.cpp:

#define EIGEN_USE_BLAS
#include <iostream>
#include "Eigen/Dense" // #include <Eigen/Dense>

int main() {
    Eigen::MatrixXd A(4, 4);
    A << 2, 1, 1, 0,
         1, 3, 1, 1,
         1, 1, 4, 1,
         0, 1, 1, 5;

    Eigen::MatrixXd A_inv = A.inverse();

    std::cout << "A:\n" << A << std::endl;
    std::cout << "\nA_inv:\n" << A_inv << std::endl;
    std::cout << "\nA*A_inv:\n" << A * A_inv << std::endl;

    return 0;
}

链接编译(只选一个):

g++ -o test openblas_test_with_eigen.cpp -lopenblas
g++ -o test openblas_test_with_eigen.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread
g++ -static -o test openblas_test_with_eigen.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread

如果编译出现报错,可能是因为编译器默认使用的是 C++98 标准,而软件包要求 C++11 或更高标准。这里考虑指定使用 C++14 标准来编译 C++ 代码,一般会解决问题(只选一个):

g++ -std=c++14 -o test openblas_test_with_eigen.cpp -lopenblas
g++ -std=c++14 -o test openblas_test_with_eigen.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread
g++ -std=c++14 -static -o test openblas_test_with_eigen.cpp ~/OpenBLAS/lib/libopenblas.a -lpthread

编译后运行:

./test

运行结果:

A:
2 1 1 0
1 3 1 1
1 1 4 1
0 1 1 5

A_inv:
  0.666667       -0.2  -0.133333  0.0666667
      -0.2       0.44      -0.04      -0.08
 -0.133333      -0.04   0.306667 -0.0533333
 0.0666667      -0.08 -0.0533333   0.226667

A*A_inv:
           1 -2.08167e-17            0  6.93889e-18
-1.38778e-16            1            0 -2.77556e-17
-1.38778e-16            0            1            0
           0            0 -5.55112e-17            1
4 次浏览

【说明:本站主要是个人的一些笔记和代码分享,内容可能会不定期修改。为了使全网显示的始终是最新版本,这里的文章未经同意请勿转载。引用请注明出处:https://www.guanjihuan.com

发表评论

您的邮箱地址不会被公开。 必填项已用 * 标注

Captcha Code