help-gsl
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Help-gsl] memory issues with gsl matrix and vector views - causes a cor


From: Ted Byers
Subject: [Help-gsl] memory issues with gsl matrix and vector views - causes a core dump/access violation
Date: Sun, 25 Mar 2012 22:32:51 -0400

The code is direct simple.

 

void regressionPCA::Reset(const std::vector<std::vector<double> >&data, 

                             const std::vector<double>& Ydata) {

  unsigned int r(0),c(0),n(0);

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  r = data.size();

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  c = data[0].size();

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  setNrows(r);

  setNcols(c);

  n = r * c;

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  B.reset(new double[c]);

  Y.reset(new double[r]); 

  U.reset(new double[n]);

  V.reset(new double[c*c]);

  S.reset(new double[c]);

  std::cout << "Flag 1" << std::endl << std::flush;

  double* Btmp = B.get();

  double* Ytmp = Y.get();

  double* Utmp = U.get();

  double* Vtmp = V.get();

  double* Stmp = S.get();

  double *bptr = Utmp;

  std::vector<std::vector<double> >::const_iterator it = data.begin(), end =
data.end();

  while (it != end) {

    bptr = std::copy(it->begin(), it->end(),bptr);

    ++it;

  }

  bptr = Ytmp;

  std::cout << "Flag 2" << std::endl << std::flush;

  std::copy(Ydata.begin(),Ydata.end(),bptr);

  std::cout << "Flag 3" << std::endl << std::flush;

  gsl_matrix_view Um = gsl_matrix_view_array(Utmp, getNrows(), getNcols());

  gsl_vector_view Ym = gsl_vector_view_array(Ytmp, getNrows());

  gsl_vector_view Bm = gsl_vector_view_array(Btmp, getNcols());

  gsl_vector_view Sm = gsl_vector_view_array(Stmp, getNcols());

  gsl_matrix_view Vm = gsl_matrix_view_array(Vtmp, getNcols(), getNcols());

  gsl_linalg_SV_decomp_jacobi(&Um.matrix,&Vm.matrix,&Sm.vector);

 
gsl_linalg_SV_solve(&Um.matrix,&Vm.matrix,&Sm.vector,&Ym.vector,&Bm.vector);

  std::cout << std::endl << std::endl << "Sv = " <<
gsl_vector_get(&Sm.vector,0) << "\t" <<  gsl_vector_get(&Sm.vector,1) <<
std::endl << std::flush;

  std::cout << std::endl << std::endl << "V = " << std::endl;

  std::cout << "\t" << gsl_matrix_get(&Vm.matrix,0,0) << "\t" <<
gsl_matrix_get(&Vm.matrix,0,1) << std::endl;

  std::cout << "\t" << gsl_matrix_get(&Vm.matrix,1,0) << "\t" <<
gsl_matrix_get(&Vm.matrix,1,1) << std::endl;

  std::cout << std::endl << std::endl << "Beta = " <<
gsl_vector_get(&Bm.vector,0) << "\t" <<  gsl_vector_get(&Bm.vector,1) <<
std::endl << std::flush;

};

 

This code executes fine the first two or three times the function it
invoked, but it never gets past the third execution without crashing.  Now,
I thought that maybe it was my use of the boost::shared_array that might be
at fault, but the following program proves that that is not the case (as it
runs to completion).

 

void regressionPCA::Reset(const std::vector<std::vector<double> >&data, 

                             const std::vector<double>& Ydata) {

  unsigned int r(0),c(0),n(0);

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  r = data.size();

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  c = data[0].size();

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  setNrows(r);

  setNcols(c);

  n = r * c;

  std::cout << "r: " << r << "\tc: " << c << "\tn: " << n << std::endl <<
std::flush;

  B.reset(new double[c]);

  Y.reset(new double[r]); 

  U.reset(new double[n]);

  V.reset(new double[c*c]);

  S.reset(new double[c]);

  std::cout << "Flag 1" << std::endl << std::flush;

  double* Btmp = B.get();

  double* Ytmp = Y.get();

  double* Utmp = U.get();

  double* Vtmp = V.get();

  double* Stmp = S.get();

  double *bptr = Utmp;

  std::vector<std::vector<double> >::const_iterator it = data.begin(), end =
data.end();

  while (it != end) {

    bptr = std::copy(it->begin(), it->end(),bptr);

    ++it;

  }

  bptr = Ytmp;

  std::cout << "Flag 2" << std::endl << std::flush;

  std::copy(Ydata.begin(),Ydata.end(),bptr);

  std::cout << "Flag 3" << std::endl << std::flush;

  gsl_matrix_view Um = gsl_matrix_view_array(Utmp, getNrows(), getNcols());

  gsl_vector_view Ym = gsl_vector_view_array(Ytmp, getNrows());

  gsl_vector_view Bm = gsl_vector_view_array(Btmp, getNcols());

  gsl_vector_view Sm = gsl_vector_view_array(Stmp, getNcols());

  gsl_matrix_view Vm = gsl_matrix_view_array(Vtmp, getNcols(), getNcols());

  gsl_linalg_SV_decomp_jacobi(&Um.matrix,&Vm.matrix,&Sm.vector);

 
gsl_linalg_SV_solve(&Um.matrix,&Vm.matrix,&Sm.vector,&Ym.vector,&Bm.vector);

  std::cout << std::endl << std::endl << "Sv = " <<
gsl_vector_get(&Sm.vector,0) << "\t" <<  gsl_vector_get(&Sm.vector,1) <<
std::endl << std::flush;

  std::cout << std::endl << std::endl << "V = " << std::endl;

  std::cout << "\t" << gsl_matrix_get(&Vm.matrix,0,0) << "\t" <<
gsl_matrix_get(&Vm.matrix,0,1) << std::endl;

  std::cout << "\t" << gsl_matrix_get(&Vm.matrix,1,0) << "\t" <<
gsl_matrix_get(&Vm.matrix,1,1) << std::endl;

  std::cout << std::endl << std::endl << "Beta = " <<
gsl_vector_get(&Bm.vector,0) << "\t" <<  gsl_vector_get(&Bm.vector,1) <<
std::endl << std::flush;

};

 

So, I went further, and if I replace the shared_array by naked pointers, and
apply operator delete[] to them at the very end of the function, the same
problem arises.  So what can possibly be wrong in this code, or is it the
case GSL is doing something odd behind the scenes that causes the core dump?

 

Cheers

 

Ted

 



reply via email to

[Prev in Thread] Current Thread [Next in Thread]