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 = 16;
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); fn compute_column_means(&self) -> Vec<T>;
27 fn multiply_with_dense(&self, dense: &DMatrix<T>, result: &mut DMatrix<T>, transpose_self: bool);
28
29 fn multiply_with_dense_centered(&self, dense: &DMatrix<T>, result: &mut DMatrix<T>, transpose_self: bool, means: &DVector<T>);
30}
31
32#[derive(Debug, Clone, PartialEq)]
41pub struct SvdRec<T: Float> {
42 pub d: usize,
43 pub u: Array2<T>,
44 pub s: Array1<T>,
45 pub vt: Array2<T>,
46 pub diagnostics: Diagnostics<T>,
47}
48
49#[derive(Debug, Clone, PartialEq)]
64pub struct Diagnostics<T: Float> {
65 pub non_zero: usize,
66 pub dimensions: usize,
67 pub iterations: usize,
68 pub transposed: bool,
69 pub lanczos_steps: usize,
70 pub ritz_values_stabilized: usize,
71 pub significant_values: usize,
72 pub singular_values: usize,
73 pub end_interval: [T; 2],
74 pub kappa: T,
75 pub random_seed: u32,
76}
77
78pub trait SvdFloat: FloatOpsTS {
79 fn eps() -> Self;
80 fn eps34() -> Self;
81 fn compare(a: Self, b: Self) -> bool;
82}
83
84impl SvdFloat for f32 {
85 fn eps() -> Self {
86 f32::EPSILON
87 }
88
89 fn eps34() -> Self {
90 f32::EPSILON.powf(0.75)
91 }
92
93 fn compare(a: Self, b: Self) -> bool {
94 (b - a).abs() < f32::EPSILON
95 }
96}
97
98impl SvdFloat for f64 {
99 fn eps() -> Self {
100 f64::EPSILON
101 }
102
103 fn eps34() -> Self {
104 f64::EPSILON.powf(0.75)
105 }
106
107 fn compare(a: Self, b: Self) -> bool {
108 (b - a).abs() < f64::EPSILON
109 }
110}