Eigen LDLT Cholesky-Zerlegung an Ort und Stelle

Question

Ich versuche, Eigen3 zu bekommen, um ein lineares System A * X = B mit einer In-Place-Cholesky-Zerlegung zu lösen. Ich kann es mir nicht leisten, irgendwelche Provisorien der Größe A auf den Stapel zu schieben, aber ich bin frei, A dabei zu zerstören. Leider

,

A.llt().solveInPlace(B); 

steht außer Frage, da A.llt() implizit eine temporäre Matrix der Größe von A auf den Stapel schiebt. Für den LLT Fall konnte ich Zugriff auf die notwendige Funktionalität erhalten wie so:

// solve A * X = B in-place for positive-definite A 
template <typename AType, typename BType> 
void AllInPlaceSolve(AType& A, BType& B) 
{ 
    typedef Eigen::internal::LLT_Traits<AType, Eigen::Upper> TraitsType; 
    TraitsType::inplace_decomposition(A); 
    TraitsType::getL(A).solveInPlace(B); 
    TraitsType::getU(A).solveInPlace(B); 
} 

Dies funktioniert gut, aber ich bin besorgt, dass:

• Meine Matrizen A könnten nur positiv semidefinit sein, in denen Fall Zersetzung ein LDLT erforderlich
• die LLT Zersetzung berechnet sqrt() unnötig für die Lösung des Systems

Ich konnte keine Möglichkeit finden, die LDLT-Funktionalität von Eigen ähnlich dem obigen Code einzubinden, da der Code sehr unterschiedlich strukturiert ist.

Meine Frage ist also: Gibt es eine Möglichkeit, Eigen3 zu verwenden, um ein lineares System zu lösen, das LDLT-Zerlegungen verwendet und keinen Platz mehr für die Diagonalmatrix D verwendet?

Answer 1

Eine Möglichkeit ist, nur einmal einen LDLT Löser zuzuweisen, und die Rechenmethode aufrufen:

LDLT<MatType> ldlt(size); 
// ... 
ldlt.compute(A); 
x = ldlt.solve(b); 

Wenn das auch nicht möglich ist, können Sie die von der LDLT Objekt gespeichert Matrix gegossen konst:

LDLT<MatType> ldlt(MatType::Identity(size,size)); 
MatType& A = const_cast<MatType&>(ldlt.matrixLDLT()); 

spielt mit A, und dann:

ldlt.compute(A); 
x = ldlt.solve(b); 

Dies ist hässlich, aber diese sollte so lange funktionieren, wie MatType Spaltenhaupt ist.