This bugzilla service is closed. All entries have been migrated to https://gitlab.com/libeigen/eigen
Bug 1381 - Segfault of expression on sparse matrices
Summary: Segfault of expression on sparse matrices
Status: RESOLVED FIXED
Alias: None
Product: Eigen
Classification: Unclassified
Component: Core - expression templates (show other bugs)
Version: 3.3 (current stable)
Hardware: All All
: Normal Crash
Assignee: Nobody
URL:
Whiteboard:
Keywords:
Depends on:
Blocks:
 
Reported: 2017-01-23 00:24 UTC by Yixuan Qiu
Modified: 2019-12-04 16:44 UTC (History)
3 users (show)



Attachments

Description Yixuan Qiu 2017-01-23 00:24:20 UTC
Hi All,

Below is a simple example that generates a segment fault on a sparse matrix expression.


=====================================================================

#include <Eigen/Core>
#include <Eigen/SparseCore>

typedef Eigen::SparseMatrix<double> SpMat;
using Eigen::ArrayXd;

int main()
{
    SpMat m(10, 10);
    ArrayXd xx = m.diagonal().array();

    return 0;
}

=====================================================================

After tracing back the call stack, I found that it was due to a return value of a local address.

In the assignment operation (https://bitbucket.org/eigen/eigen/src/aeba495019d15b6cb90f9ae2455e22d56f542fee/Eigen/src/Core/functors/AssignmentFunctors.h?at=3.3&fileviewer=file-view-default#AssignmentFunctors.h-24), the signature of the source data is "const SrcScalar& b". However, the sparse matrix expression can only return a local "Scalar" (https://bitbucket.org/eigen/eigen/src/aeba495019d15b6cb90f9ae2455e22d56f542fee/Eigen/src/SparseCore/SparseCompressedBase.h?at=3.3&fileviewer=file-view-default#SparseCompressedBase.h-300). Therefore, a local address will be returned (https://bitbucket.org/eigen/eigen/src/aeba495019d15b6cb90f9ae2455e22d56f542fee/Eigen/src/Core/CoreEvaluators.h?at=3.3&fileviewer=file-view-default#CoreEvaluators.h-1340) to the assignment operator that results a reading data error.


Best,
Yixuan
Comment 1 Gael Guennebaud 2017-01-25 15:19:03 UTC
On your side you can omit the ".array()".
Comment 2 Gael Guennebaud 2017-01-25 16:41:19 UTC
Fixed:

https://bitbucket.org/eigen/eigen/commits/f5b44011f514/
https://bitbucket.org/eigen/eigen/commits/44d70389b325/ (3.3)
Summary:     Bug 1381: fix sparse.diagonal() used as a rvalue.
The problem was that is "sparse" is not const, then sparse.diagonal() must have the
LValueBit flag meaning that sparse.diagonal().coeff(i) must returns a const reference,
const Scalar&. However, sparse::coeff() cannot returns a reference for a non-existing
zero coefficient. The trick is to return a reference to a local member of
evaluator<SparseMatrix>.
Comment 3 Yixuan Qiu 2017-01-25 18:08:09 UTC
Thank you!

The background is that I'm preparing for an update of RcppEigen (R wrapper of Eigen), and I checked the compatibility of Eigen 3.3 with the dependent R packages. One package failed due to this bug, so I created a minimal example to mimic its behavior.
Comment 4 Nobody 2019-12-04 16:44:53 UTC
-- GitLab Migration Automatic Message --

This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity.

You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/1381.

Note You need to log in before you can comment on or make changes to this bug.