Is there any memory issue with the way I handled this?

Question

I am writing a wrapper to Eigen QR for my personal use and I am wondering if there is any memory leak or undocumented behavior in my implementation especially in the function void get_QR(double* A, int m, int n, double*& Q, double*& R). The answer is as expected. This is related to my previous question here.

using std::cout;
using std::endl;
using namespace Eigen;

/*!
 Obtains the QR decomposition as A=QR, where all the matrices are in Eigen MatrixXd format.
 */
void get_QR(MatrixXd A, MatrixXd& Q, MatrixXd& R) {

        int m           =       A.rows();
        int n           =       A.cols();
        int minmn       =       min(m,n);

        //      A_E     =       Q_E*R_E.
        HouseholderQR<MatrixXd> qr(A);
        Q  =   qr.householderQ()*(MatrixXd::Identity(m, minmn));
        R  =   qr.matrixQR().block(0, 0, minmn, n).triangularView<Upper>();
}

/*!
 Obtains the QR decomposition as A=QR, where all the matrices are in double format.
 */
void get_QR(double* A, int m, int n, double*& Q, double*& R) {
        MatrixXd Q_E, R_E;
        int minmn       =       min(m,n);

        //      Maps the double to MatrixXd.
        Map<MatrixXd> A_E(A, m, n);

        get_QR(A_E, Q_E, R_E);

        Q       =       (double*)realloc(Q_E.data(), m*minmn*sizeof(double));
        R       =       (double*)realloc(R_E.data(), minmn*n*sizeof(double));
}


int main(int argc, char* argv[]) {
        srand(time(NULL));

        int m   =       atoi(argv[1]);
        int n   =       atoi(argv[2]);

//      Check the double version.
        int minmn       =       min(m,n);
        double* A       =       (double*)malloc(m*n*sizeof(double));
        double* Q       =       (double*)malloc(m*minmn*sizeof(double));
        double* R       =       (double*)malloc(minmn*n*sizeof(double));

        double RANDMAX  =       double(RAND_MAX);
        //      Initialize A as a random matrix.
        for (int index=0; index<m*n; ++index) {
                A[index]        =       rand()/RANDMAX;
        }

        get_QR(A, m, n, Q, R);
        std::cout << Q[0] << std::endl;

//      Check the MatrixXd version.
        Map<MatrixXd> A_E(A, m, n);
        MatrixXd Q_E, R_E;
        get_QR(A_E, Q_E, R_E);

        cout << Q[0] << endl;
        cout << Q_E(0,0) << endl;

        free(A);
        free(Q);
        free(R);
}

For instance, I get the output as

-0.360995
-0.360995
-0.360995

Answer 1

You can replace your calls to malloc and free with a vector. Here is sample code of how you would use vector with existing code that takes a pointer to double.

#include <vector>
//...
void foo(double* d, int size)
{
   d[0] = 10;
}

int main()
{
   std::vector<double> d(10);  // 10 doubles
   foo(&d[0], d.size());
}

The foo() function is an "old C" function that expects a pointer to double to denote a double array. Note that to do this magic, you pass the address of the first element.

So your original code could look like this, given the above:

#include <vector>
#include <algorithm>

typedef std::vector<double> VectDouble;

void get_QR(VectDouble& A, int m, int n, VectDouble& Q, VectDouble& R) 
{
    MatrixXd Q_E, R_E;
    int minmn       =       min(m,n);

    //      Maps the double to MatrixXd.
    Map<MatrixXd> A_E(&A[0], m, n);

    get_QR(A_E, Q_E, R_E);

    Q.resize(m * minmn);
    R.resize(minmn *n, 0);
}

int main(int argc, char* argv[]) 
{
    srand(time(NULL));

    int m   =       atoi(argv[1]);
    int n   =       atoi(argv[2]);

    int minmn       =       min(m,n);
    VectDouble A(m*n);
    VectDouble Q(m*minmn);
    VectDouble R(minmn*n);

    double RANDMAX  =       double(RAND_MAX);
    //      Initialize A as a random matrix.
    std::fill(A.begin(), A.end(), rand()/RANDMAX);

    get_QR(A, m, n, Q, R);
    std::cout << Q[0] << std::endl;

    Map<MatrixXd> A_E(&A[0], m, n);
    MatrixXd Q_E, R_E;
    get_QR(A_E, Q_E, R_E);

    cout << Q[0] << endl;
    cout << Q_E(0,0) << endl;
}

Answer 2

realloc is not guaranteed to be compatible with operator new[] and operator delete[], only with malloc and free.

I would just replace all your double* pointers with std::vector<double> objects.