I'm trying to write an interior point solver for linear optimization problems. I've run into three issues I haven't been able to answer by reading the docs or forum topics. (Warning: I'm pretty new to both C++ and Eigen, so I could be missing simple things.)

(1) Is it a good idea to use (C++11) smart pointers with MatrixXf* objects? I want to use them to manage some MatrixXf data members, but I'm concerned about the cost of the pointer de-reference in formulas, or if using the pointers might interfere with the expression templates.

(2) When I use the (column pivot) QR decomposition for a matrix A, this seems to return an object that contains all the factors. Is this correct? I'm also unsure what the simplest way to get those factors without allocating more memory is. It seems like the type of the returned QR object is, for a MatrixXf type, Eigen::ColPivHouseholderQR<Eigen::MatrixXf>, and I get Q via householderQ(), R via matrixQR().triangularView<Eigen::Upper>(), and P via colsPermutation().

Is there a simpler way to access the factors?

(3) The interior point solver is iterative. Some variables are replaced from iteration to iteration, but there's no reason the memory allocated for them couldn't be re-used. (The vector/matrix sizes are only determined at run time instead of compile-time, and will be large in any case. So I think I have to use MatrixXf.) Since that memory is already allocated, is it possible to re-use it from iteration to iteration? So, for example, if b is a vector with three blocks, will the expression template magic make it so

Code: Select all

v.segment(0,0+n) = b-A*x;


will not create any temporaries?

(4) The compiler wasn't giving me errors when I accessed clearly out of bounds elements, such as

Code: Select all

MatrixXf m(3,3);
m(3,3) = 5;


It also didn't give an error when I multiplied a permutation matrix and a MatrixXf of non-compatible sizes. I didn't compile with any flags to disable debug mode, so I was surprised that it failed at runtime rather than compile time. What can I do to make the compiler give an error in these cases?

Thanks!

(1) There should not be issues, but check twice that smart pointers are really worth it in your case. Is the lifetime of objects that complex to justify smart pointers?

(2) What do you mean by simpler ?

(3) Instead of:

v.segment(0,0+n) = b-A*x;

write:

v.segment(0,0+n) = b:
v.segment(0,0+n) .noalias() -= A*x;

to get rid of a temporary.

(4) such errors can only be detected at runtime, and they will be if you don't compile with -DNDEBUG or -DEIGEN_NO_DEBUG