diff --git a/src/linalg/cholesky.rs b/src/linalg/cholesky.rs index 0362f809..7623e247 100644 --- a/src/linalg/cholesky.rs +++ b/src/linalg/cholesky.rs @@ -6,8 +6,9 @@ use alga::general::ComplexField; use crate::allocator::Allocator; use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, SquareMatrix}; use crate::constraint::{SameNumberOfRows, ShapeConstraint}; -use crate::dimension::{Dim, DimSub, Dynamic, U1}; +use crate::dimension::{Dim, DimAdd, DimSum, DimSub, Dynamic, U1}; use crate::storage::{Storage, StorageMut}; +use crate::base::allocator::Reallocator; /// The Cholesky decomposition of a symmetric-definite-positive matrix. #[cfg_attr(feature = "serde-serialize", derive(Serialize, Deserialize))] @@ -188,6 +189,34 @@ where } } } + + /// Updates the decomposition such that we get the decomposition of a matrix with the given column `c` in the `j`th position. + /// Since the matrix is square, an identical row will be added in the `j`th row. + pub fn insert_column( + self, + j: usize, + c: &Matrix, + ) -> Cholesky> + where + D: DimAdd, + DefaultAllocator: Reallocator> + Reallocator, DimSum, DimSum>, + S2: Storage, + ShapeConstraint: SameNumberOfRows>, + { + let n = c.nrows(); + assert_eq!( + n, + self.chol.nrows() + 1, + "The new column must have the size of the factored matrix plus one." + ); + assert!(j < n, "j needs to be within the bound of the new matrix."); + // TODO what is the fastest way to produce the new matrix ? + let chol= self.chol.insert_column(j, N::zero()).insert_row(j, N::zero()); + + // TODO see https://en.wikipedia.org/wiki/Cholesky_decomposition#Updating_the_decomposition + unimplemented!(); + Cholesky { chol } + } } impl, S: Storage> SquareMatrix