pub struct Cholesky<T: SimdComplexField, D: Dim>where
DefaultAllocator: Allocator<D, D>,{ /* private fields */ }
Expand description
The Cholesky decomposition of a symmetric-definite-positive matrix.
Implementations§
source§impl<T: SimdComplexField, D: Dim> Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
impl<T: SimdComplexField, D: Dim> Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
sourcepub fn new_unchecked(matrix: OMatrix<T, D, D>) -> Self
pub fn new_unchecked(matrix: OMatrix<T, D, D>) -> Self
Computes the Cholesky decomposition of matrix
without checking that the matrix is definite-positive.
If the input matrix is not definite-positive, the decomposition may contain trash values (Inf, NaN, etc.)
sourcepub fn pack_dirty(matrix: OMatrix<T, D, D>) -> Self
pub fn pack_dirty(matrix: OMatrix<T, D, D>) -> Self
Uses the given matrix as-is without any checks or modifications as the Cholesky decomposition.
It is up to the user to ensure all invariants hold.
sourcepub fn unpack(self) -> OMatrix<T, D, D>
pub fn unpack(self) -> OMatrix<T, D, D>
Retrieves the lower-triangular factor of the Cholesky decomposition with its strictly upper-triangular part filled with zeros.
sourcepub fn unpack_dirty(self) -> OMatrix<T, D, D>
pub fn unpack_dirty(self) -> OMatrix<T, D, D>
Retrieves the lower-triangular factor of the Cholesky decomposition, without zeroing-out its strict upper-triangular part.
The values of the strict upper-triangular part are garbage and should be ignored by further computations.
sourcepub fn l(&self) -> OMatrix<T, D, D>
pub fn l(&self) -> OMatrix<T, D, D>
Retrieves the lower-triangular factor of the Cholesky decomposition with its strictly uppen-triangular part filled with zeros.
sourcepub fn l_dirty(&self) -> &OMatrix<T, D, D>
pub fn l_dirty(&self) -> &OMatrix<T, D, D>
Retrieves the lower-triangular factor of the Cholesky decomposition, without zeroing-out its strict upper-triangular part.
This is an allocation-less version of self.l()
. The values of the strict upper-triangular
part are garbage and should be ignored by further computations.
sourcepub fn solve_mut<R2: Dim, C2: Dim, S2>(&self, b: &mut Matrix<T, R2, C2, S2>)
pub fn solve_mut<R2: Dim, C2: Dim, S2>(&self, b: &mut Matrix<T, R2, C2, S2>)
Solves the system self * x = b
where self
is the decomposed matrix and x
the unknown.
The result is stored on b
.
sourcepub fn solve<R2: Dim, C2: Dim, S2>(
&self,
b: &Matrix<T, R2, C2, S2>
) -> OMatrix<T, R2, C2>where
S2: Storage<T, R2, C2>,
DefaultAllocator: Allocator<R2, C2>,
ShapeConstraint: SameNumberOfRows<R2, D>,
pub fn solve<R2: Dim, C2: Dim, S2>(
&self,
b: &Matrix<T, R2, C2, S2>
) -> OMatrix<T, R2, C2>where
S2: Storage<T, R2, C2>,
DefaultAllocator: Allocator<R2, C2>,
ShapeConstraint: SameNumberOfRows<R2, D>,
Returns the solution of the system self * x = b
where self
is the decomposed matrix and
x
the unknown.
sourcepub fn determinant(&self) -> T::SimdRealField
pub fn determinant(&self) -> T::SimdRealField
Computes the determinant of the decomposed matrix.
sourcepub fn ln_determinant(&self) -> T::SimdRealField
pub fn ln_determinant(&self) -> T::SimdRealField
Computes the natural logarithm of determinant of the decomposed matrix.
This method is more robust than .determinant()
to very small or very
large determinants since it returns the natural logarithm of the
determinant rather than the determinant itself.
source§impl<T: ComplexField, D: Dim> Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
impl<T: ComplexField, D: Dim> Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
sourcepub fn new(matrix: OMatrix<T, D, D>) -> Option<Self>
pub fn new(matrix: OMatrix<T, D, D>) -> Option<Self>
Attempts to compute the Cholesky decomposition of matrix
.
Returns None
if the input matrix is not definite-positive. The input matrix is assumed
to be symmetric and only the lower-triangular part is read.
sourcepub fn new_with_substitute(
matrix: OMatrix<T, D, D>,
substitute: T
) -> Option<Self>
pub fn new_with_substitute( matrix: OMatrix<T, D, D>, substitute: T ) -> Option<Self>
Attempts to approximate the Cholesky decomposition of matrix
by
replacing non-positive values on the diagonals during the decomposition
with the given substitute
.
try_sqrt
will be applied to the substitute
when it has to be used.
If your input matrix results only in positive values on the diagonals
during the decomposition, substitute
is unused and the result is just
the same as if you used new
.
This method allows to compensate for matrices with very small or even
negative values due to numerical errors but necessarily results in only
an approximation: it is basically a hack. If you don’t specifically need
Cholesky, it may be better to consider alternatives like the
LU
decomposition/factorization.
sourcepub fn rank_one_update<R2: Dim, S2>(
&mut self,
x: &Vector<T, R2, S2>,
sigma: T::RealField
)where
S2: Storage<T, R2, U1>,
DefaultAllocator: Allocator<R2, U1>,
ShapeConstraint: SameNumberOfRows<R2, D>,
pub fn rank_one_update<R2: Dim, S2>(
&mut self,
x: &Vector<T, R2, S2>,
sigma: T::RealField
)where
S2: Storage<T, R2, U1>,
DefaultAllocator: Allocator<R2, U1>,
ShapeConstraint: SameNumberOfRows<R2, D>,
Given the Cholesky decomposition of a matrix M
, a scalar sigma
and a vector v
,
performs a rank one update such that we end up with the decomposition of M + sigma * (v * v.adjoint())
.
sourcepub fn insert_column<R2, S2>(
&self,
j: usize,
col: Vector<T, R2, S2>
) -> Cholesky<T, DimSum<D, U1>>
pub fn insert_column<R2, S2>( &self, j: usize, col: Vector<T, R2, S2> ) -> Cholesky<T, DimSum<D, U1>>
Updates the decomposition such that we get the decomposition of a matrix with the given column col
in the j
th position.
Since the matrix is square, an identical row will be added in the j
th row.
Trait Implementations§
source§impl<T: Clone + SimdComplexField, D: Clone + Dim> Clone for Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
impl<T: Clone + SimdComplexField, D: Clone + Dim> Clone for Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
source§impl<T: Debug + SimdComplexField, D: Debug + Dim> Debug for Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
impl<T: Debug + SimdComplexField, D: Debug + Dim> Debug for Cholesky<T, D>where
DefaultAllocator: Allocator<D, D>,
impl<T: SimdComplexField, D: Dim> Copy for Cholesky<T, D>
Auto Trait Implementations§
impl<T, D> !Freeze for Cholesky<T, D>
impl<T, D> !RefUnwindSafe for Cholesky<T, D>
impl<T, D> !Send for Cholesky<T, D>
impl<T, D> !Sync for Cholesky<T, D>
impl<T, D> !Unpin for Cholesky<T, D>
impl<T, D> !UnwindSafe for Cholesky<T, D>
Blanket Implementations§
source§impl<T> BorrowMut<T> for Twhere
T: ?Sized,
impl<T> BorrowMut<T> for Twhere
T: ?Sized,
source§fn borrow_mut(&mut self) -> &mut T
fn borrow_mut(&mut self) -> &mut T
source§impl<SS, SP> SupersetOf<SS> for SPwhere
SS: SubsetOf<SP>,
impl<SS, SP> SupersetOf<SS> for SPwhere
SS: SubsetOf<SP>,
source§fn to_subset(&self) -> Option<SS>
fn to_subset(&self) -> Option<SS>
self
from the equivalent element of its
superset. Read moresource§fn is_in_subset(&self) -> bool
fn is_in_subset(&self) -> bool
self
is actually part of its subset T
(and can be converted to it).source§fn to_subset_unchecked(&self) -> SS
fn to_subset_unchecked(&self) -> SS
self.to_subset
but without any property checks. Always succeeds.source§fn from_subset(element: &SS) -> SP
fn from_subset(element: &SS) -> SP
self
to the equivalent element of its superset.