I've successfully narrowed down my confusion to a complete yet simple example below. Hopefully, this will help illustrate my specific confusion so you can increase my understanding of what is really going on inside the computer. It appears my current code is off writing to unallocated memory addresses in the vector and matrix. I've started with an easy example then built on it below. Please see my specific questions and exception inline in the comments. Thanks for your help and clarification. The power and speed of C++ come from the ability to do just about anything as I adjust what I think I told the computer to do with what it is really doing.
// multiply matrix and vector allocating size in constructor
void test1(int n)
{
boost::numeric::ublas::compressed_matrix<float> M(n,n);
boost::numeric::ublas::vector<float> v1(n);
boost::numeric::ublas::vector<float> v2(n);
// fill in some values on diagonal for testing...
for (int i = 0; i < n; i++)
{
M(i,i) = i + 1;
v1(i) = i + 1;
}
noalias(v2) = prod(M,v1);
}
// multiply matrix and vector dynamically allocating size
void test2(int n)
{
boost::numeric::ublas::compressed_matrix<float> *M;
boost::numeric::ublas::vector<float> *v1;
boost::numeric::ublas::vector<float> *v2;
M = new boost::numeric::ublas::compressed_matrix<float>(n,n);
v1 = new boost::numeric::ublas::vector<float>(n);
v2 = new boost::numeric::ublas::vector<float>(n);
for (int i = 0; i < n; i++)
{
(*M)(i,i) = i + 1;
(*v1)(i) = i + 1;
}
noalias(*v2) = prod(*M,*v1);
delete M;
delete v1;
delete v2;
}
/*
multiply matrix and vector attempting to resize existing objects
What does this error message mean?
Check failed in file ...boost_1_54_0/boost/numeric/ublas/matrix_sparse.hpp at line 2807:
!preserve
libc++abi.dylib: terminate called throwing an exception
*/
void test3(int n)
{
boost::numeric::ublas::compressed_matrix<float> M; // size zero
boost::numeric::ublas::vector<float> v1; // size zero
boost::numeric::ublas::vector<float> v2; // size zero
// cannot call reserve--do not know how many non zero elements will be?
// what is the correct way to set size of existing matrix and vector?
M.resize(n,n); // *** Crash?!? see exception in comment above
v1.resize(n,false); // do not preserve data
v2.resize(n,false);
for (int i = 0; i < n; i++)
{
M(i,i) = i + 1;
v1(i) = i + 1;
}
noalias(v2) = prod(M,v1);
}
// multiply vector of matrix by vector of vectors
void test4(int n1, int n2)
{
// vector[i] faster while vector.at(i) checks for array out of bounds exception
std::vector<boost::numeric::ublas::compressed_matrix<float>> theMatrix; // size zero
std::vector<boost::numeric::ublas::vector<float>> v1; // size zero
std::vector<boost::numeric::ublas::vector<float>> v2; // size zero
// all matrix and vectors same size for this test
// need to have each one a different size at runtime
// is there an easier/faster/more efficient way than dynamically shown in test5?
boost::numeric::ublas::compressed_matrix<float> M(n2,n2); // copied into vector
boost::numeric::ublas::vector<float> v(n2); // copied into vector
theMatrix.assign(n1,M); // allocate and populate container
v1.assign(n1,v); // allocate and populate container
v2.assign(n1,v); // allocate and populate container
for (int i = 0; i < n1; i++)
{
for (int j = 0; j < n2; j++)
{
theMatrix.at(i)(j,j) = j + 1;
v1.at(i)(j) = j + 1;
}
noalias(v2.at(i)) = prod(theMatrix.at(i),v1.at(i));
}
}
// multiply dynamically allocated vector of matrix by vector of vectors
// extra credit: write similar test with smart pointers from std or boost?
void test5(int n1, int n2)
{
std::vector<boost::numeric::ublas::compressed_matrix<float>*> theMatrix;
std::vector<boost::numeric::ublas::vector<float>*> v1;
std::vector<boost::numeric::ublas::vector<float>*> v2;
boost::numeric::ublas::compressed_matrix<float> *M;
boost::numeric::ublas::vector<float> *v;
for (int i = 0; i < n1; i++)
{
M = new boost::numeric::ublas::compressed_matrix<float>(n2,n2);
theMatrix.push_back(M);
v = new boost::numeric::ublas::vector<float>(n2);
v1.push_back(v);
v = new boost::numeric::ublas::vector<float>(n2);
v2.push_back(v); // need to allocate space for result--this object replaced by result from prod?
for (int j = 0; j < n2; j++)
{
(*theMatrix.at(i))(j,j) = j + 1;
(*v1.at(i))(j) = j + 1;
}
// what happens to the memory of v2 here?
noalias(*v2.at(i)) = prod(*theMatrix.at(i),*v1.at(i));
}
// while the vectors clean up themselves, the pointers do not free automatically
std::vector<boost::numeric::ublas::compressed_matrix<float>*>::iterator mit;
for (mit = theMatrix.begin(); mit != theMatrix.end(); ++mit)
{
delete *mit;
}
std::vector<boost::numeric::ublas::vector<float>*>::iterator it1;
for (it1 = v1.begin(); it1 != v1.end(); ++it1)
{
delete *it1;
}
std::vector<boost::numeric::ublas::vector<float>*>::iterator it2;
for (it2 = v2.begin(); it2 != v2.end(); ++it2)
{
delete *it1;
}
}