犰狳 eig_sym 中的线程安全问题

Question

我在 eig_sym 中对犰狳的特征分解有疑问。当我尝试并行计算多组特征值和特征向量时，特征向量有时是

不正交
未规范化
甚至没有相关矩阵的特征向量。

如果每次只有一个计算运行，这个问题就消失了（所以这似乎是一些线程安全问题）。很快，当两个计算运行并行时，问题又出现了。奇怪的是，特征值似乎在每种情况下都是正确的。

//compile with: g++ -std=c++11 -pthread -larmadillo -o evecs armadillo_evecs.cpp
#include <iostream>
#include <armadillo>
#include <assert.h>
#include <future>
#include <vector>
using namespace std;

void testcalc() {
    // set up random symmetric matrix
    arma::mat r = arma::randu<arma::mat> (100, 100);
    r  = r.t() * r;

    arma::vec eval;
    arma::mat evec;
    // calculate eigenvalues and -vectors
    assert(arma::eig_sym(eval, evec, r));
    arma::mat test = evec.t() * evec;

    // Check whether eigenvectors are orthogonal, (i. e. matrix 'test' is diagonal)
    assert(arma::norm(test - arma::diagmat(test)) < 1.0e-10);
}


int main() {
    // start 100 eigenvalue (+vector) calculations
    vector<future<void>> fus;
    for (size_t i = 0; i < 100; i++) {
        // try parallel evaluation ... fails sometimes 
        fus.push_back(async(launch::async, &testcalc));

        // try sequential evaluation ... works fine
        // future<void> f = async(launch::async, &testcalc);
        // f.get(); // Wait until calculation has finished, before starting new one
    }

    // wait until calculations have finished
    for(auto it = fus.begin(); it != fus.end(); it++) {
        it->get();
    }

    return 0;
}

所以在断言上面的代码中

assert(arma::norm(test - arma::diagmat(test)) < 1.0e-10);

有时会失败。这可能是底层库的问题（我读到 lapack 有一些线程安全问题）？我真的不知道，从哪里开始寻找。

Answer 1

而不是"rolling your own parallelization"，使用底层库已经提供的并行化更容易也更安全。

因此，不要使用参考 BLAS 和 LAPACK，而是使用像 OpenBLAS or Intel MKL. See the Armadillo FAQ 这样的多线程版本以获得更多详细信息。

犰狳 eig_sym 中的线程安全问题

Thread-safety issue in armadillo's eig_sym

c++

thread-safety

armadillo

eigenvalue

eigenvector