single_svdlib/
utils.rs

1use rayon::iter::ParallelIterator;
2use nalgebra_sparse::na::{DMatrix, DVector};
3use ndarray::{Array1, Array2};
4use num_traits::{Float, Zero};
5use rayon::prelude::{IntoParallelIterator, IndexedParallelIterator};
6use single_utilities::traits::FloatOpsTS;
7use std::fmt::Debug;
8
9pub fn determine_chunk_size(nrows: usize) -> usize {
10    let num_threads = rayon::current_num_threads();
11
12    let min_rows_per_thread = 64;
13    let desired_chunks_per_thread = 4;
14
15    let target_total_chunks = num_threads * desired_chunks_per_thread;
16    let chunk_size = nrows.div_ceil(target_total_chunks);
17
18    chunk_size.max(min_rows_per_thread)
19}
20
21pub trait SMat<T: Float> {
22    fn nrows(&self) -> usize;
23    fn ncols(&self) -> usize;
24    fn nnz(&self) -> usize;
25    fn svd_opa(&self, x: &[T], y: &mut [T], transposed: bool); // y = A*x
26    fn compute_column_means(&self) -> Vec<T>;
27    fn multiply_with_dense(&self, dense: &DMatrix<T>, result: &mut DMatrix<T>, transpose_self: bool);
28    fn multiply_with_dense_centered(&self, dense: &DMatrix<T>, result: &mut DMatrix<T>, transpose_self: bool, means: &DVector<T>);
29
30    fn multiply_transposed_by_dense(&self, q: &DMatrix<T>, result: &mut DMatrix<T>);
31    fn multiply_transposed_by_dense_centered(&self, q: &DMatrix<T>, result: &mut DMatrix<T>, means: &DVector<T>);
32}
33
34/// Singular Value Decomposition Components
35///
36/// # Fields
37/// - d:  Dimensionality (rank), the number of rows of both `ut`, `vt` and the length of `s`
38/// - ut: Transpose of left singular vectors, the vectors are the rows of `ut`
39/// - s:  Singular values (length `d`)
40/// - vt: Transpose of right singular vectors, the vectors are the rows of `vt`
41/// - diagnostics: Computational diagnostics
42#[derive(Debug, Clone, PartialEq)]
43pub struct SvdRec<T: Float> {
44    pub d: usize,
45    pub u: Array2<T>,
46    pub s: Array1<T>,
47    pub vt: Array2<T>,
48    pub diagnostics: Diagnostics<T>,
49}
50
51/// Computational Diagnostics
52///
53/// # Fields
54/// - non_zero:  Number of non-zeros in the matrix
55/// - dimensions: Number of dimensions attempted (bounded by matrix shape)
56/// - iterations: Number of iterations attempted (bounded by dimensions and matrix shape)
57/// - transposed:  True if the matrix was transposed internally
58/// - lanczos_steps: Number of Lanczos steps performed
59/// - ritz_values_stabilized: Number of ritz values
60/// - significant_values: Number of significant values discovered
61/// - singular_values: Number of singular values returned
62/// - end_interval: left, right end of interval containing unwanted eigenvalues
63/// - kappa: relative accuracy of ritz values acceptable as eigenvalues
64/// - random_seed: Random seed provided or the seed generated
65#[derive(Debug, Clone, PartialEq)]
66pub struct Diagnostics<T: Float> {
67    pub non_zero: usize,
68    pub dimensions: usize,
69    pub iterations: usize,
70    pub transposed: bool,
71    pub lanczos_steps: usize,
72    pub ritz_values_stabilized: usize,
73    pub significant_values: usize,
74    pub singular_values: usize,
75    pub end_interval: [T; 2],
76    pub kappa: T,
77    pub random_seed: u32,
78}
79
80pub trait SvdFloat: FloatOpsTS {
81    fn eps() -> Self;
82    fn eps34() -> Self;
83    fn compare(a: Self, b: Self) -> bool;
84}
85
86impl SvdFloat for f32 {
87    fn eps() -> Self {
88        f32::EPSILON
89    }
90
91    fn eps34() -> Self {
92        f32::EPSILON.powf(0.75)
93    }
94
95    fn compare(a: Self, b: Self) -> bool {
96        (b - a).abs() < f32::EPSILON
97    }
98}
99
100impl SvdFloat for f64 {
101    fn eps() -> Self {
102        f64::EPSILON
103    }
104
105    fn eps34() -> Self {
106        f64::EPSILON.powf(0.75)
107    }
108
109    fn compare(a: Self, b: Self) -> bool {
110        (b - a).abs() < f64::EPSILON
111    }
112}
single_svdlib/utils.rs

single_svdlib/
utils.rs
