I am wanting to find the nearest positive definite matrix to some X matrix.
I have seen this answer here: https://scicomp.stackexchange.com/questions/30631/how-to-find-the-nearest-a-near-positive-definite-from-a-given-matrix
But I am having trouble translating that to an answer in Eigen. I have tried the following:
Matrix<double, 9, 9> X, D, DPLUS, Z, Y;
X.setRandom();
Y = 0.5 * (X + X.transpose().eval());
EigenSolver<Matrix<double, 9, 9>> es;
es.compute(Y);
D = es.eigenvalues().asDiagonal();
DPlus = D.cwiseMax(0);
Z = es.eigenvectors() * DPLUS * es.eigenvectors().transpose().eval();
This is giving me errors about a complex problem, but is this inline with what the linked answer is suggesting?
The problems that you're not using the right eigenvalue solver.
Since the matrix Y is symmetric by construction, you need to use SelfAdjointEigenSolver.
This symmetry guarantees that your real matrix is diagonalizable and has real eigenvalues and eigenvectors, so you don't have to deal with complex numbers.
The general EigenSolver doesn't make this assumption and thus isn't going to work as you intend in this case. You could get the same result with it (you'll get complex eigen values which are all going to be in fact reals, but you'll need to create vectors of complexes etc. Also it is much less efficient).
This will give you the eigen values and the eigen vectors which you can use to construct your orthogonal basis matrix Q.
const MatrixXd Y = 0.5 * (X + X.transpose());
const SelfAdjointEigenSolver<MatrixXd> solver(Y);
const VectorXd D = solver.eigenvalues();
const MatrixXd Q = solver.eigenvectors();
const VectorXd Dplus = D.cwiseMax(0);
const MatrixXd Z = Q * Dplus.asDiagonal() * Q.transpose();
Of course if you like one-liners you can always shorten it to :
const SelfAdjointEigenSolver<MatrixXd> solver(0.5 * (X + X.transpose()));
const MatrixXd Z = solver.eigenvectors() * solver.eigenvalues().cwiseMax(0).asDiagonal() * solver.eigenvectors().transpose();