srenblad/nalgebra
1
0
Fork 0
forked from M-Labs/nalgebra
Code Pull Requests Activity

tests pass, needs cleanup

Browse Source
This commit is contained in:
Nestor Demeure 2019-11-03 18:48:04 +01:00
parent e6467c2a3e
commit 46687b7cde
2 changed files with 9 additions and 6 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
src/linalg/cholesky.rs
View File

@ -216,18 +216,21 @@ where
let mut chol= self.chol.clone().insert_column(j, N::zero()).insert_row(j, N::zero()); let mut chol= self.chol.clone().insert_column(j, N::zero()).insert_row(j, N::zero());
// update the jth row // update the jth row
let top_left_corner = chol.slice_range(..j-1, ..j-1); let top_left_corner = chol.slice_range(..j, ..j);
let colj_minus = col.rows_range(..j-1); let colj_minus = col.rows_range(..j);
let rowj = top_left_corner.solve_lower_triangular(&colj_minus).unwrap().adjoint(); // TODO both the row and its adjoint seem to be usefull let rowj = top_left_corner.solve_lower_triangular(&colj_minus).unwrap().adjoint(); // TODO both the row and its adjoint seem to be usefull
chol.slice_range_mut(j, ..j-1).copy_from(&rowj); chol.slice_range_mut(j, ..j).copy_from(&rowj);
// TODO
//println!("dotc:{} norm2:{}", rowj.dotc(&rowj), rowj.norm_squared());
// update the center element // update the center element
let center_element = N::sqrt(col[j] + rowj.dot(&rowj.adjoint())); // TODO is there a better way to multiply a vector by its adjoint ? norm_squared ? let center_element = N::sqrt(col[j] - rowj.dotc(&rowj) );
chol[(j,j)] = center_element; chol[(j,j)] = center_element;
// update the jth column // update the jth column
let colj_plus = col.rows_range(j+1..); let colj_plus = col.rows_range(j+1..);
let bottom_left_corner = chol.slice_range(j+1.., ..j-1); let bottom_left_corner = chol.slice_range(j+1.., ..j);
let colj = (colj_plus - bottom_left_corner*rowj.adjoint()) / center_element; // TODO that can probably be done with a single optimized operation let colj = (colj_plus - bottom_left_corner*rowj.adjoint()) / center_element; // TODO that can probably be done with a single optimized operation
chol.slice_range_mut(j+1.., j).copy_from(&colj); chol.slice_range_mut(j+1.., j).copy_from(&colj);

2
tests/linalg/cholesky.rs
View File

@ -100,7 +100,7 @@ macro_rules! gen_tests(
} }
fn cholesky_insert_column(n: usize) -> bool { fn cholesky_insert_column(n: usize) -> bool {
let n = n.max(1).min(5); let n = n.max(1).min(50);
let j = random::<usize>() % n; let j = random::<usize>() % n;
let m_updated = RandomSDP::new(Dynamic::new(n), || random::<$scalar>().0).unwrap(); let m_updated = RandomSDP::new(Dynamic::new(n), || random::<$scalar>().0).unwrap();