[eigen] Memory fault when using a function returning a SelfAdjointView

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]


In statistical computing I will frequently want to calculate X'X from
an m by n model matrix X where m >= n.  I use the rankUpdate method
for the SelfAdjointView<MatrixXd, Lower> class to calculate the
symmetric result.

I can do this inline but if I try to calculate the value in a function
returning the SelfAdjointView, rather than a MatrixXd, I encounter a
bad_alloc error under Ubuntu 11.10 (64-bit).  I enclose an example
program.

The reason that I want to return the SelfAdjointView is because I want
to create different types of decompositions, such as LLT, LDLT and
SelfAdjointEigenSolver from the result.  If I return the value as a
matrix I think I need to apply the selfadjointView method before
creating the decomposition, which seems redundant.
#include <iostream>
#include <eigen3/Eigen/Dense>

using namespace std;
using namespace Eigen;

SelfAdjointView<MatrixXd, Lower> XtX(const MatrixXd& X) {
    return MatrixXd(X.cols(), X.cols()).setZero().selfadjointView<Lower>().rankUpdate(X.adjoint());
}

MatrixXd  XtXmat(const MatrixXd& X) {
    return MatrixXd(MatrixXd(X.cols(), X.cols()).setZero().
		    selfadjointView<Lower>().rankUpdate(X.adjoint()));
}

int main() {
    MatrixXd X(4, 2);

    X <<
	1, 1,
	1, 2,
	1, 3,
	1, 4;

    // inline calculation of XtX is okay
    MatrixXd  xtx(MatrixXd(X.cols(), X.cols()).setZero().
		  selfadjointView<Lower>().rankUpdate(X.adjoint()));
    cout << "XtX(inline):\n" << xtx << endl;

    // returning XtX as a MatrixXd is okay
    cout << "XtX(mat):\n" << XtXmat(X) << endl;

    // but not returning a SelfAdjointView

    MatrixXd  XtXres(XtX(X));
    cout << "XtX(function):\n" << XtXres << endl;
}


Mail converted by MHonArc 2.6.19+ http://listengine.tuxfamily.org/