blu/update.rs
1// Copyright (C) 2016-2018 ERGO-Code
2// Copyright (C) 2022-2023 Richard Lincoln
3
4use crate::lu;
5use crate::lu::LU;
6use crate::Status;
7
8/// Update the factorization to replace one column of the factorized matrix.
9/// A call to [`update()`] must be preceded by calls to
10/// [`solve_for_update()`](crate::solve_for_update()) to provide the column
11/// to be inserted and the index of the column to be replaced.
12///
13/// The column to be inserted is defined as the right-hand side in the last call
14/// to [`solve_for_update()`](crate::solve_for_update()) in which the forward
15/// system was solved.
16///
17/// The index of the column to be replaced is defined by the unit vector in the
18/// last call to [`solve_for_update()`](crate::solve_for_update()) in which the
19/// transposed system was solved.
20///
21/// ## Arguments
22///
23/// `xtbl` is an optional argument to monitor numerical stability. `xtbl` can be
24/// either of:
25///
26/// - element `j0` of the solution to the forward system computed by
27/// [`solve_for_update()`](crate::solve_for_update()), where `j0` is the column
28/// to be replaced;
29/// - the dot product of the incoming column and the solution to the transposed
30/// system computed by [`solve_for_update()`](crate::solve_for_update()).
31///
32/// In either case [`LU::pivot_error`] has a defined value. If monitoring
33/// stability is not desired, `xtbl` can be any value.
34///
35/// ## Returns
36///
37/// [`Status::ErrorInvalidCall`] if the factorization is invalid or the update was not
38/// prepared by two calls to [`solve_for_update()`](crate::solve_for_update()).
39///
40/// [`Status::Reallocate`] for insufficient memory in `w_i`,`w_x`. The number of additional
41/// elements required is given by [`LU::addmem_w`]. The user must reallocate
42/// `w_i`,`w_x`. It is recommended to reallocate for the requested number of additional
43/// elements plus some extra space for further updates (e.g. 0.5 times the current array length).
44/// The new array length must be provided in [`LU::w_mem`]. [`update()`] will start
45/// from scratch in the next call.
46///
47/// [`Status::ErrorSingularUpdate`] if the updated factorization would be (numerically) singular.
48/// No update has been computed and the old factorization is still valid.
49pub fn update(lu: &mut LU, xtbl: f64) -> Status {
50 if lu.nupdate.is_none() || lu.ftran_for_update.is_none() || lu.btran_for_update.is_none() {
51 Status::ErrorInvalidCall
52 } else {
53 lu::update(lu, xtbl)
54 }
55}