CODEWITHSUNDEEP
c++thread-safetyarmadilloeigenvalueeigenvector
okruz
okruz

Reputation: 95

Thread-safety issue in armadillo's eig_sym

I have a problem with the eigen decomposition in armadillo via eig_sym. When I try to calculate multiple sets of eigenvalues and eigenvectors in parallel, from time to time the eigenvectors are

This problem vanishes, if only one calculation is run at each time (so this seems to be some thread-safety issue). As soon, as two calculations run in parallel, the problem crops up again. Strangely, the eigenvalues seem to be correct in every case.

//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;
}

So in the code above the assertion 

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

fails sometimes. Might this be a problem of the underlying libraries (I read lapack has had some thread-safety issues)? I don't really know, where to start looking.

Upvotes: 1

Views: 488

Answers (1)

hbrerkere
hbrerkere

Reputation: 1615

Instead of "rolling your own parallelization", it's easier and safer to use the parallelization already provided by an underlying library.

So instead of using reference BLAS and LAPACK, use a multi-threaded version like OpenBLAS or Intel MKL. See the Armadillo FAQ for more details.

Upvotes: 1

Related Questions