new_identity and housholder matrix available under na::

This commit is contained in:
Vincent Barrielle 2014-05-12 21:54:59 +02:00
parent 3dc0e27fd2
commit 73c6610048
3 changed files with 16 additions and 11 deletions

View File

@ -11,11 +11,11 @@ use std::cmp::min;
/// * `dim` - the dimension of the space the resulting matrix operates in /// * `dim` - the dimension of the space the resulting matrix operates in
/// * `start` - the starting dimension of the subspace of the reflexion /// * `start` - the starting dimension of the subspace of the reflexion
/// * `vec` - the vector defining the reflection. /// * `vec` - the vector defining the reflection.
fn householder_matrix<N: Float, pub fn householder_matrix<N: Float,
Mat: Eye + Indexable<(uint, uint), N>, M: Eye + Indexable<(uint, uint), N>,
V: Indexable<uint, N>> V: Indexable<uint, N>>
(dim: uint, start: uint, vec: V) -> Mat { (dim: uint, start: uint, vec: V) -> M {
let mut qk : Mat = Eye::new_identity(dim); let mut qk : M = Eye::new_identity(dim);
let stop = start + vec.shape(); let stop = start + vec.shape();
assert!(stop <= dim); assert!(stop <= dim);
for j in range(start, stop) { for j in range(start, stop) {
@ -35,12 +35,12 @@ fn householder_matrix<N: Float,
/// * `m` - matrix to decompose /// * `m` - matrix to decompose
pub fn decomp_qr<N: Float, pub fn decomp_qr<N: Float,
V: Indexable<uint, N> + Norm<N>, V: Indexable<uint, N> + Norm<N>,
Mat: Clone + Eye + ColSlice<V> + Transpose M: Clone + Eye + ColSlice<V> + Transpose
+ Indexable<(uint, uint), N> + Mul<Mat, Mat>> + Indexable<(uint, uint), N> + Mul<M, M>>
(m: &Mat) -> (Mat, Mat) { (m: &M) -> (M, M) {
let (rows, cols) = m.shape(); let (rows, cols) = m.shape();
assert!(rows >= cols); assert!(rows >= cols);
let mut q : Mat = Eye::new_identity(rows); let mut q : M = Eye::new_identity(rows);
let mut r = m.clone(); let mut r = m.clone();
let iterations = min(rows - 1, cols); let iterations = min(rows - 1, cols);
@ -59,7 +59,7 @@ pub fn decomp_qr<N: Float,
v.unsafe_set(0, x - alpha); v.unsafe_set(0, x - alpha);
} }
let _ = v.normalize(); let _ = v.normalize();
let qk: Mat = householder_matrix(rows, 0, v); let qk: M = householder_matrix(rows, 0, v);
r = qk * r; r = qk * r;
q = q * Transpose::transpose_cpy(&qk); q = q * Transpose::transpose_cpy(&qk);
} }

View File

@ -1,4 +1,4 @@
pub use self::decompositions::decomp_qr; pub use self::decompositions::{decomp_qr, householder_matrix};
mod decompositions; mod decompositions;

View File

@ -55,7 +55,8 @@ pub use structs::{
}; };
pub use linalg::{ pub use linalg::{
decomp_qr decomp_qr,
householder_matrix
}; };
/// Traits to work around the language limitations related to operator overloading. /// Traits to work around the language limitations related to operator overloading.
@ -721,6 +722,10 @@ pub fn mean<N, M: Mean<N>>(observations: &M) -> N {
// //
// //
/// Construct the identity matrix for a given dimension
#[inline(always)]
pub fn new_identity<M: Eye>(dim: uint) -> M { Eye::new_identity(dim) }
/* /*
* Basis * Basis
*/ */