scirs2_sparse/
lib.rs

1#![allow(clippy::manual_div_ceil)]
2#![allow(clippy::needless_return)]
3#![allow(clippy::manual_ok_err)]
4#![allow(clippy::needless_range_loop)]
5#![allow(clippy::while_let_loop)]
6#![allow(clippy::vec_init_then_push)]
7#![allow(clippy::should_implement_trait)]
8#![allow(clippy::only_used_in_recursion)]
9#![allow(clippy::manual_slice_fill)]
10#![allow(dead_code)]
11//! # SciRS2 Sparse - Sparse Matrix Operations
12//!
13//! **scirs2-sparse** provides comprehensive sparse matrix formats and operations modeled after SciPy's
14//! `sparse` module, offering CSR, CSC, COO, DOK, LIL, DIA, BSR formats with efficient algorithms
15//! for large-scale sparse linear algebra, eigenvalue problems, and graph operations.
16//!
17//! ## 🎯 Key Features
18//!
19//! - **SciPy Compatibility**: Drop-in replacement for `scipy.sparse` classes
20//! - **Multiple Formats**: CSR, CSC, COO, DOK, LIL, DIA, BSR with easy conversion
21//! - **Efficient Operations**: Sparse matrix-vector/matrix multiplication
22//! - **Linear Solvers**: Direct (LU, Cholesky) and iterative (CG, GMRES) solvers
23//! - **Eigenvalue Solvers**: ARPACK-based sparse eigenvalue computation
24//! - **Array API**: Modern NumPy-compatible array interface (recommended)
25//!
26//! ## 📦 Module Overview
27//!
28//! | SciRS2 Format | SciPy Equivalent | Description |
29//! |---------------|------------------|-------------|
30//! | `CsrArray` | `scipy.sparse.csr_array` | Compressed Sparse Row (efficient row slicing) |
31//! | `CscArray` | `scipy.sparse.csc_array` | Compressed Sparse Column (efficient column slicing) |
32//! | `CooArray` | `scipy.sparse.coo_array` | Coordinate format (efficient construction) |
33//! | `DokArray` | `scipy.sparse.dok_array` | Dictionary of Keys (efficient element access) |
34//! | `LilArray` | `scipy.sparse.lil_array` | List of Lists (efficient incremental construction) |
35//! | `DiaArray` | `scipy.sparse.dia_array` | Diagonal format (efficient banded matrices) |
36//! | `BsrArray` | `scipy.sparse.bsr_array` | Block Sparse Row (efficient block operations) |
37//!
38//! ## 🚀 Quick Start
39//!
40//! ```toml
41//! [dependencies]
42//! scirs2-sparse = "0.4.3"
43//! ```
44//!
45//! ```rust
46//! use scirs2_sparse::csr_array::CsrArray;
47//!
48//! // Create sparse matrix from triplets (row, col, value)
49//! let rows = vec![0, 0, 1, 2, 2];
50//! let cols = vec![0, 2, 2, 0, 1];
51//! let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
52//! let sparse = CsrArray::from_triplets(&rows, &cols, &data, (3, 3), false).expect("Operation failed");
53//! ```
54//!
55//! ## 🔒 Version: 0.4.3 (March 27, 2026)
56//!
57//! ## Matrix vs. Array API
58//!
59//! This module provides both a matrix-based API and an array-based API,
60//! following SciPy's transition to a more NumPy-compatible array interface.
61//!
62//! When using the array interface (e.g., `CsrArray`), please note that:
63//!
64//! - `*` performs element-wise multiplication, not matrix multiplication
65//! - Use `dot()` method for matrix multiplication
66//! - Operations like `sum` produce arrays, not matrices
67//! - Array-style slicing operations return scalars, 1D, or 2D arrays
68//!
69//! For new code, we recommend using the array interface, which is more consistent
70//! with the rest of the numerical ecosystem.
71//!
72//! ## Examples
73//!
74//! ### Matrix API (Legacy)
75//!
76//! ```
77//! use scirs2_sparse::csr::CsrMatrix;
78//!
79//! // Create a sparse matrix in CSR format
80//! let rows = vec![0, 0, 1, 2, 2];
81//! let cols = vec![0, 2, 2, 0, 1];
82//! let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
83//! let shape = (3, 3);
84//!
85//! let matrix = CsrMatrix::new(data, rows, cols, shape).expect("Operation failed");
86//! ```
87//!
88//! ### Array API (Recommended)
89//!
90//! ```
91//! use scirs2_sparse::csr_array::CsrArray;
92//!
93//! // Create a sparse array in CSR format
94//! let rows = vec![0, 0, 1, 2, 2];
95//! let cols = vec![0, 2, 2, 0, 1];
96//! let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
97//! let shape = (3, 3);
98//!
99//! // From triplets (COO-like construction)
100//! let array = CsrArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
101//!
102//! // Or directly from CSR components
103//! // let array = CsrArray::new(...);
104//! ```
105
106// Export error types
107pub mod error;
108pub use error::{SparseError, SparseResult};
109
110// Base trait for sparse arrays
111pub mod sparray;
112pub use sparray::{is_sparse, SparseArray, SparseSum};
113
114// Trait for symmetric sparse arrays
115pub mod sym_sparray;
116pub use sym_sparray::SymSparseArray;
117
118// No spatial module in sparse
119
120// Array API (recommended)
121pub mod csr_array;
122pub use csr_array::CsrArray;
123
124pub mod csc_array;
125pub use csc_array::CscArray;
126
127pub mod coo_array;
128pub use coo_array::CooArray;
129
130pub mod dok_array;
131pub use dok_array::DokArray;
132
133pub mod lil_array;
134pub use lil_array::LilArray;
135
136pub mod dia_array;
137pub use dia_array::DiaArray;
138
139pub mod bsr_array;
140pub use bsr_array::BsrArray;
141
142pub mod banded_array;
143pub use banded_array::BandedArray;
144
145// Symmetric array formats
146pub mod sym_csr;
147pub use sym_csr::{SymCsrArray, SymCsrMatrix};
148
149pub mod sym_coo;
150pub use sym_coo::{SymCooArray, SymCooMatrix};
151
152// Legacy matrix formats
153pub mod csr;
154pub use csr::CsrMatrix;
155
156pub mod csc;
157pub use csc::CscMatrix;
158
159pub mod coo;
160pub use coo::CooMatrix;
161
162pub mod dok;
163pub use dok::DokMatrix;
164
165pub mod lil;
166pub use lil::LilMatrix;
167
168pub mod dia;
169pub use dia::DiaMatrix;
170
171pub mod bsr;
172pub use bsr::BsrMatrix;
173
174pub mod banded;
175pub use banded::BandedMatrix;
176
177// Utility functions
178pub mod utils;
179
180// Linear algebra with sparse matrices
181pub mod linalg;
182// Re-export the main functions from the reorganized linalg module
183pub use linalg::{
184    // Functions from solvers
185    add,
186    // Functions from iterative
187    bicg,
188    bicgstab,
189    cg,
190    cholesky_decomposition,
191    // Enhanced operators
192    convolution_operator,
193    diag_matrix,
194    eigs,
195    eigsh,
196    enhanced_add,
197    enhanced_diagonal,
198    enhanced_scale,
199    enhanced_subtract,
200    expm,
201    // Functions from matfuncs
202    expm_multiply,
203    eye,
204    finite_difference_operator,
205    // GCROT solver
206    gcrot,
207    gmres,
208    incomplete_cholesky,
209    incomplete_lu,
210    inv,
211    // IRAM eigenvalue solver (v0.3.0)
212    iram,
213    iram_shift_invert,
214    lanczos,
215    // Decomposition functions
216    lu_decomposition,
217    matmul,
218    matrix_power,
219    multiply,
220    norm,
221    onenormest,
222    // Eigenvalue functions
223    power_iteration,
224    qr_decomposition,
225    // Specialized solvers (v0.2.0)
226    solve_arrow_matrix,
227    solve_banded_system,
228    solve_block_2x2,
229    solve_kronecker_system,
230    solve_saddle_point,
231    sparse_direct_solve,
232    sparse_lstsq,
233    spsolve,
234    svd_truncated,
235    // SVD functions
236    svds,
237    // TFQMR solver
238    tfqmr,
239    // IRAM (Arnoldi) eigenvalue solver (v0.3.0)
240    ArnoldiConfig,
241    ArpackOptions,
242    // Interfaces
243    AsLinearOperator,
244    // Types from iterative
245    BiCGOptions,
246    BiCGSTABOptions,
247    BiCGSTABResult,
248    // Enhanced operator types
249    BoundaryCondition,
250    CGOptions,
251    CGSOptions,
252    CGSResult,
253    CholeskyResult,
254    ConvolutionMode,
255    ConvolutionOperator,
256    // Operator types
257    DiagonalOperator,
258    EigenResult,
259    EigenvalueMethod,
260    EnhancedDiagonalOperator,
261    EnhancedDifferenceOperator,
262    EnhancedOperatorOptions,
263    EnhancedScaledOperator,
264    EnhancedSumOperator,
265    FiniteDifferenceOperator,
266    GCROTOptions,
267    GCROTResult,
268    GMRESOptions,
269    ICOptions,
270    // Preconditioners
271    ILU0Preconditioner,
272    ILUOptions,
273    IdentityOperator,
274    IterationResult,
275    JacobiPreconditioner,
276    // Decomposition types
277    LUResult,
278    LanczosOptions,
279    LinearOperator,
280    // Eigenvalue types
281    PowerIterationOptions,
282    QRResult,
283    SSORPreconditioner,
284    // SVD types
285    SVDOptions,
286    SVDResult,
287    ScaledIdentityOperator,
288    TFQMROptions,
289    TFQMRResult,
290};
291
292// Format conversions
293pub mod convert;
294
295// Construction utilities
296pub mod construct;
297pub mod construct_sym;
298
299// Combining arrays
300pub mod combine;
301pub use combine::{block_diag, bmat, hstack, kron, kronsum, tril, triu, vstack};
302
303// Index dtype handling utilities
304pub mod index_dtype;
305pub use index_dtype::{can_cast_safely, get_index_dtype, safely_cast_index_arrays};
306
307// Optimized operations for symmetric sparse formats
308pub mod sym_ops;
309
310// Tensor-based sparse operations (v0.2.0)
311pub mod tensor_sparse;
312
313// Enhanced BSR format with flat storage and Block LU (v0.3.0)
314pub mod bsr_enhanced;
315pub use bsr_enhanced::{block_lu, block_lu_solve, BlockLUResult, EnhancedBsr};
316
317// Enhanced DIA format with banded operations (v0.3.0)
318pub mod dia_enhanced;
319pub use dia_enhanced::{banded_solve, dia_tridiagonal_solve, tridiagonal_solve, EnhancedDia};
320
321// Compressed Sparse Fiber (CSF) tensor format (v0.3.0)
322pub mod csf_tensor;
323pub use csf_tensor::CsfTensor;
324
325// Advanced sparse formats (v0.4.0): SELL-C-sigma, CSR5, CSF (standalone), Polynomial Preconditioners
326pub mod formats;
327pub use formats::csf::CsfTensor as CsfTensorStandalone;
328pub use formats::csr5::Csr5Matrix;
329pub use formats::poly_precond::{
330    ChebyshevPreconditioner, NeumannPreconditioner as NeumannPreconditionerPoly, PolyPrecondConfig,
331};
332pub use formats::sell::SellMatrix;
333
334// Sparse matrix utility functions (v0.3.0)
335pub mod sparse_functions;
336pub use sparse_functions::{
337    sparse_block_diag, sparse_diag_matrix, sparse_diags, sparse_eye, sparse_eye_rect,
338    sparse_hstack, sparse_kron, sparse_kronsum, sparse_random, sparse_vstack,
339};
340pub use sym_ops::{
341    sym_coo_matvec, sym_csr_matvec, sym_csr_quadratic_form, sym_csr_rank1_update, sym_csr_trace,
342};
343
344// Tensor operations (v0.2.0)
345pub use tensor_sparse::{khatri_rao_product, CPDecomposition, SparseTensor, TuckerDecomposition};
346
347// GPU-accelerated operations
348pub mod gpu;
349pub mod gpu_kernel_execution;
350pub mod gpu_ops;
351pub mod gpu_preconditioner;
352pub mod gpu_spmv_implementation;
353pub use gpu_kernel_execution::{
354    calculate_adaptive_workgroup_size, execute_spmv_kernel, execute_symmetric_spmv_kernel,
355    execute_triangular_solve_kernel, GpuKernelConfig, GpuMemoryManager as GpuKernelMemoryManager,
356    GpuPerformanceProfiler, MemoryStrategy,
357};
358pub use gpu_ops::{
359    gpu_sparse_matvec, gpu_sym_sparse_matvec, AdvancedGpuOps, GpuKernelScheduler, GpuMemoryManager,
360    GpuOptions, GpuProfiler, OptimizedGpuOps,
361};
362pub use gpu_spmv_implementation::GpuSpMV;
363
364// Memory-efficient algorithms and patterns
365pub mod memory_efficient;
366pub use memory_efficient::{
367    streaming_sparse_matvec, CacheAwareOps, MemoryPool, MemoryTracker, OutOfCoreProcessor,
368};
369
370// SIMD-accelerated operations
371pub mod simd_ops;
372pub use simd_ops::{
373    simd_csr_matvec, simd_sparse_elementwise, simd_sparse_linear_combination, simd_sparse_matmul,
374    simd_sparse_norm, simd_sparse_scale, simd_sparse_transpose, ElementwiseOp, SimdOptions,
375};
376
377// Parallel vector operations for iterative solvers
378pub mod parallel_vector_ops;
379pub use parallel_vector_ops::{
380    advanced_sparse_matvec_csr, parallel_axpy, parallel_dot, parallel_linear_combination,
381    parallel_norm2, parallel_sparse_matvec_csr, parallel_vector_add, parallel_vector_copy,
382    parallel_vector_scale, parallel_vector_sub, ParallelVectorOptions,
383};
384
385// Enhanced iterative solvers with preconditioners and sparse utilities (v0.3.0)
386pub mod iterative_solvers;
387pub use iterative_solvers::{
388    // Solvers
389    bicgstab as enhanced_bicgstab,
390    cg as enhanced_cg,
391    chebyshev,
392    // Sparse utility functions
393    estimate_spectral_radius,
394    gmres as enhanced_gmres,
395    sparse_diagonal,
396    sparse_norm,
397    sparse_trace,
398    // Preconditioners (Array1-based interface)
399    ILU0Preconditioner as EnhancedILU0Preconditioner,
400    // Configuration and result types
401    IterativeSolverConfig,
402    JacobiPreconditioner as EnhancedJacobiPreconditioner,
403    NormType,
404    Preconditioner,
405    SSORPreconditioner as EnhancedSSORPreconditioner,
406    SolverResult,
407};
408
409// LOBPCG eigensolver (v0.3.0)
410pub mod lobpcg;
411pub use lobpcg::{
412    lobpcg as lobpcg_eigensolver, lobpcg_generalised, EigenTarget, LobpcgConfig, LobpcgResult,
413};
414
415// Advanced Krylov subspace eigensolvers (v0.3.0)
416pub mod krylov;
417pub use krylov::{
418    iram as krylov_iram, thick_restart_lanczos, IramConfig, KrylovEigenResult,
419    ThickRestartLanczosConfig, WhichEigenvalues,
420};
421
422// Sparse matrix utilities (v0.3.0)
423pub mod sparse_utils;
424pub use sparse_utils::{
425    condest_1norm, permute_matrix, reverse_cuthill_mckee, sparse_add, sparse_extract_diagonal,
426    sparse_identity, sparse_kronecker, sparse_matrix_norm, sparse_matrix_trace, sparse_scale,
427    sparse_sub, sparse_transpose, spgemm, RcmResult, SparseNorm,
428};
429
430// Incomplete factorization preconditioners (v0.3.0)
431pub mod incomplete_factorizations;
432pub use incomplete_factorizations::{Ic0, Ilu0 as Ilu0Enhanced, IluK, Ilut, IlutConfig, Milu};
433
434// Sparse direct solvers (v0.3.0)
435pub mod direct_solver;
436pub use direct_solver::{
437    amd_ordering, elimination_tree, inverse_perm, nested_dissection_ordering,
438    sparse_cholesky_solve, sparse_lu_solve, symbolic_cholesky, SparseCholResult,
439    SparseCholeskySolver, SparseLuResult, SparseLuSolver, SparseSolver, SymbolicAnalysis,
440};
441
442// Sparse QR factorization (v0.3.0)
443pub mod sparse_qr;
444pub use sparse_qr::{
445    apply_q, apply_qt, extract_q_dense, numerical_rank, sparse_least_squares,
446    sparse_qr as sparse_qr_factorize, SparseLeastSquaresResult, SparseQrConfig, SparseQrResult,
447};
448
449// Sparse matrix reordering algorithms (v0.4.0)
450pub mod reorder;
451pub use reorder::{
452    amd, amd_simple, apply_permutation as reorder_apply_permutation, apply_permutation_csr_array,
453    bandwidth as reorder_bandwidth, cuthill_mckee, cuthill_mckee_full, distance2_color,
454    dsatur_color, greedy_color, nested_dissection, nested_dissection_full,
455    nested_dissection_with_config, profile as reorder_profile,
456    reverse_cuthill_mckee as reorder_rcm, reverse_cuthill_mckee_full, verify_coloring,
457    verify_distance2_coloring, AdjacencyGraph, AmdResult, ColoringResult, CuthillMcKeeResult,
458    GreedyOrdering, NestedDissectionConfig, NestedDissectionResult,
459};
460
461// Unified sparse eigenvalue interface (v0.3.0)
462pub mod sparse_eigen;
463pub use sparse_eigen::{
464    cayley_transform_matvec, check_eigenpairs, compute_residuals, shift_invert_eig, sparse_eig,
465    sparse_eigs, sparse_eigsh, EigenMethod, EigenvalueTarget, SparseEigenConfig, SparseEigenResult,
466    SpectralTransform,
467};
468
469// Quantum-inspired sparse matrix operations (Advanced mode)
470pub mod quantum_inspired_sparse;
471pub use quantum_inspired_sparse::{
472    QuantumProcessorStats, QuantumSparseConfig, QuantumSparseProcessor, QuantumStrategy,
473};
474
475// Neural-adaptive sparse matrix operations (Advanced mode)
476pub mod neural_adaptive_sparse;
477pub use neural_adaptive_sparse::{
478    NeuralAdaptiveConfig, NeuralAdaptiveSparseProcessor, NeuralProcessorStats, OptimizationStrategy,
479};
480
481// Quantum-Neural hybrid optimization (Advanced mode)
482pub mod quantum_neural_hybrid;
483pub use quantum_neural_hybrid::{
484    HybridStrategy, QuantumNeuralConfig, QuantumNeuralHybridProcessor, QuantumNeuralHybridStats,
485};
486
487// Adaptive memory compression for advanced-large sparse matrices (Advanced mode)
488pub mod adaptive_memory_compression;
489pub use adaptive_memory_compression::{
490    AdaptiveCompressionConfig, AdaptiveMemoryCompressor, CompressedMatrix, CompressionAlgorithm,
491    MemoryStats,
492};
493
494// Real-time performance monitoring and adaptation (Advanced mode)
495pub mod realtime_performance_monitor;
496pub use realtime_performance_monitor::{
497    Alert, AlertSeverity, PerformanceMonitorConfig, PerformanceSample, ProcessorType,
498    RealTimePerformanceMonitor,
499};
500
501// Cholesky update/downdate
502pub mod cholesky_update;
503// Distributed sparse operations (SpMV halo exchange, dist AMG)
504pub mod distributed;
505// Learned multigrid smoother
506pub mod learned_smoother;
507// Low-rank sparse update
508pub mod low_rank_update;
509// ML-based preconditioner
510pub mod ml_preconditioner;
511// Parallel AMG
512pub mod parallel_amg;
513
514// Compressed sparse graph algorithms
515pub mod csgraph;
516pub use csgraph::{
517    all_pairs_shortest_path,
518    bellman_ford_single_source,
519    // Centrality measures (v0.2.0)
520    betweenness_centrality,
521    bfs_distances,
522    // Traversal algorithms
523    breadth_first_search,
524    closeness_centrality,
525    // Community detection (v0.2.0)
526    community_detection,
527    compute_laplacianmatrix,
528    connected_components,
529    degree_matrix,
530    depth_first_search,
531    dijkstra_single_source,
532    // Max flow algorithms (v0.2.0)
533    dinic,
534    edmonds_karp,
535    eigenvector_centrality,
536    floyd_warshall,
537    ford_fulkerson,
538    has_path,
539    is_connected,
540    is_laplacian,
541    is_spanning_tree,
542    // Minimum spanning trees
543    kruskal_mst,
544    label_propagation,
545    // Laplacian matrices
546    laplacian,
547    largest_component,
548    louvain_communities,
549    min_cut,
550    minimum_spanning_tree,
551    modularity,
552    num_edges,
553    num_vertices,
554    pagerank,
555    prim_mst,
556    reachable_vertices,
557    reconstruct_path,
558    // Graph algorithms
559    shortest_path,
560    // Shortest path algorithms
561    single_source_shortest_path,
562    spanning_tree_weight,
563    strongly_connected_components,
564    to_adjacency_list,
565    topological_sort,
566    traversegraph,
567    // Connected components
568    undirected_connected_components,
569    // Graph utilities
570    validate_graph,
571    weakly_connected_components,
572    LaplacianType,
573    MSTAlgorithm,
574    // Max flow types (v0.2.0)
575    MaxFlowResult,
576    // Enums and types
577    ShortestPathMethod,
578    TraversalOrder,
579};
580
581// Graph Laplacian and spectral graph operations
582pub mod graph_laplacian;
583pub use graph_laplacian::{
584    fiedler_vector, graph_laplacian, normalized_laplacian, random_walk_laplacian,
585};
586
587// Re-export warnings from scipy for compatibility
588pub struct SparseEfficiencyWarning;
589pub struct SparseWarning;
590
591/// Check if an object is a sparse array
592#[allow(dead_code)]
593pub fn is_sparse_array<T>(obj: &dyn SparseArray<T>) -> bool
594where
595    T: scirs2_core::SparseElement + std::ops::Div<Output = T> + PartialOrd + 'static,
596{
597    sparray::is_sparse(obj)
598}
599
600/// Check if an object is a symmetric sparse array
601#[allow(dead_code)]
602pub fn is_sym_sparse_array<T>(obj: &dyn SymSparseArray<T>) -> bool
603where
604    T: scirs2_core::SparseElement
605        + std::ops::Div<Output = T>
606        + scirs2_core::Float
607        + PartialOrd
608        + 'static,
609{
610    obj.is_symmetric()
611}
612
613/// Check if an object is a sparse matrix (legacy API)
614#[allow(dead_code)]
615pub fn is_sparse_matrix(obj: &dyn std::any::Any) -> bool {
616    obj.is::<CsrMatrix<f64>>()
617        || obj.is::<CscMatrix<f64>>()
618        || obj.is::<CooMatrix<f64>>()
619        || obj.is::<DokMatrix<f64>>()
620        || obj.is::<LilMatrix<f64>>()
621        || obj.is::<DiaMatrix<f64>>()
622        || obj.is::<BsrMatrix<f64>>()
623        || obj.is::<SymCsrMatrix<f64>>()
624        || obj.is::<SymCooMatrix<f64>>()
625        || obj.is::<CsrMatrix<f32>>()
626        || obj.is::<CscMatrix<f32>>()
627        || obj.is::<CooMatrix<f32>>()
628        || obj.is::<DokMatrix<f32>>()
629        || obj.is::<LilMatrix<f32>>()
630        || obj.is::<DiaMatrix<f32>>()
631        || obj.is::<BsrMatrix<f32>>()
632        || obj.is::<SymCsrMatrix<f32>>()
633        || obj.is::<SymCooMatrix<f32>>()
634}
635
636#[cfg(test)]
637mod tests {
638    use super::*;
639    use approx::assert_relative_eq;
640
641    #[test]
642    fn test_csr_array() {
643        let rows = vec![0, 0, 1, 2, 2];
644        let cols = vec![0, 2, 2, 0, 1];
645        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
646        let shape = (3, 3);
647
648        let array =
649            CsrArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
650
651        assert_eq!(array.shape(), (3, 3));
652        assert_eq!(array.nnz(), 5);
653        assert!(is_sparse_array(&array));
654    }
655
656    #[test]
657    fn test_coo_array() {
658        let rows = vec![0, 0, 1, 2, 2];
659        let cols = vec![0, 2, 2, 0, 1];
660        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
661        let shape = (3, 3);
662
663        let array =
664            CooArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
665
666        assert_eq!(array.shape(), (3, 3));
667        assert_eq!(array.nnz(), 5);
668        assert!(is_sparse_array(&array));
669    }
670
671    #[test]
672    fn test_dok_array() {
673        let rows = vec![0, 0, 1, 2, 2];
674        let cols = vec![0, 2, 2, 0, 1];
675        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
676        let shape = (3, 3);
677
678        let array = DokArray::from_triplets(&rows, &cols, &data, shape).expect("Operation failed");
679
680        assert_eq!(array.shape(), (3, 3));
681        assert_eq!(array.nnz(), 5);
682        assert!(is_sparse_array(&array));
683
684        // Test setting and getting values
685        let mut array = DokArray::<f64>::new((2, 2));
686        array.set(0, 0, 1.0).expect("Operation failed");
687        array.set(1, 1, 2.0).expect("Operation failed");
688
689        assert_eq!(array.get(0, 0), 1.0);
690        assert_eq!(array.get(0, 1), 0.0);
691        assert_eq!(array.get(1, 1), 2.0);
692
693        // Test removing zeros
694        array.set(0, 0, 0.0).expect("Operation failed");
695        assert_eq!(array.nnz(), 1);
696    }
697
698    #[test]
699    fn test_lil_array() {
700        let rows = vec![0, 0, 1, 2, 2];
701        let cols = vec![0, 2, 2, 0, 1];
702        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
703        let shape = (3, 3);
704
705        let array = LilArray::from_triplets(&rows, &cols, &data, shape).expect("Operation failed");
706
707        assert_eq!(array.shape(), (3, 3));
708        assert_eq!(array.nnz(), 5);
709        assert!(is_sparse_array(&array));
710
711        // Test setting and getting values
712        let mut array = LilArray::<f64>::new((2, 2));
713        array.set(0, 0, 1.0).expect("Operation failed");
714        array.set(1, 1, 2.0).expect("Operation failed");
715
716        assert_eq!(array.get(0, 0), 1.0);
717        assert_eq!(array.get(0, 1), 0.0);
718        assert_eq!(array.get(1, 1), 2.0);
719
720        // Test sorted indices
721        assert!(array.has_sorted_indices());
722
723        // Test removing zeros
724        array.set(0, 0, 0.0).expect("Operation failed");
725        assert_eq!(array.nnz(), 1);
726    }
727
728    #[test]
729    fn test_dia_array() {
730        use scirs2_core::ndarray::Array1;
731
732        // Create a 3x3 diagonal matrix with main diagonal + upper diagonal
733        let data = vec![
734            Array1::from_vec(vec![1.0, 2.0, 3.0]), // Main diagonal
735            Array1::from_vec(vec![4.0, 5.0, 0.0]), // Upper diagonal
736        ];
737        let offsets = vec![0, 1]; // Main diagonal and k=1
738        let shape = (3, 3);
739
740        let array = DiaArray::new(data, offsets, shape).expect("Operation failed");
741
742        assert_eq!(array.shape(), (3, 3));
743        assert_eq!(array.nnz(), 5); // 3 on main diagonal, 2 on upper diagonal
744        assert!(is_sparse_array(&array));
745
746        // Test values
747        assert_eq!(array.get(0, 0), 1.0);
748        assert_eq!(array.get(1, 1), 2.0);
749        assert_eq!(array.get(2, 2), 3.0);
750        assert_eq!(array.get(0, 1), 4.0);
751        assert_eq!(array.get(1, 2), 5.0);
752        assert_eq!(array.get(0, 2), 0.0);
753
754        // Test from_triplets
755        let rows = vec![0, 0, 1, 1, 2];
756        let cols = vec![0, 1, 1, 2, 2];
757        let data_vec = vec![1.0, 4.0, 2.0, 5.0, 3.0];
758
759        let array2 =
760            DiaArray::from_triplets(&rows, &cols, &data_vec, shape).expect("Operation failed");
761
762        // Should have same values
763        assert_eq!(array2.get(0, 0), 1.0);
764        assert_eq!(array2.get(1, 1), 2.0);
765        assert_eq!(array2.get(2, 2), 3.0);
766        assert_eq!(array2.get(0, 1), 4.0);
767        assert_eq!(array2.get(1, 2), 5.0);
768
769        // Test conversion to other formats
770        let csr = array.to_csr().expect("Operation failed");
771        assert_eq!(csr.nnz(), 5);
772        assert_eq!(csr.get(0, 0), 1.0);
773        assert_eq!(csr.get(0, 1), 4.0);
774    }
775
776    #[test]
777    fn test_format_conversions() {
778        let rows = vec![0, 0, 1, 2, 2];
779        let cols = vec![0, 2, 1, 0, 2];
780        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
781        let shape = (3, 3);
782
783        // Create a COO array
784        let coo =
785            CooArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
786
787        // Convert to CSR
788        let csr = coo.to_csr().expect("Operation failed");
789
790        // Check values are preserved
791        let coo_dense = coo.to_array();
792        let csr_dense = csr.to_array();
793
794        for i in 0..shape.0 {
795            for j in 0..shape.1 {
796                assert_relative_eq!(coo_dense[[i, j]], csr_dense[[i, j]]);
797            }
798        }
799    }
800
801    #[test]
802    fn test_dot_product() {
803        let rows = vec![0, 0, 1, 2, 2];
804        let cols = vec![0, 2, 1, 0, 2];
805        let data = vec![1.0, 2.0, 3.0, 4.0, 5.0];
806        let shape = (3, 3);
807
808        // Create arrays in different formats
809        let coo =
810            CooArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
811        let csr =
812            CsrArray::from_triplets(&rows, &cols, &data, shape, false).expect("Operation failed");
813
814        // Compute dot product (matrix multiplication)
815        let coo_result = coo.dot(&coo).expect("Operation failed");
816        let csr_result = csr.dot(&csr).expect("Operation failed");
817
818        // Check results match
819        let coo_dense = coo_result.to_array();
820        let csr_dense = csr_result.to_array();
821
822        for i in 0..shape.0 {
823            for j in 0..shape.1 {
824                assert_relative_eq!(coo_dense[[i, j]], csr_dense[[i, j]], epsilon = 1e-10);
825            }
826        }
827    }
828
829    #[test]
830    fn test_sym_csr_array() {
831        // Create a symmetric matrix
832        let data = vec![2.0, 1.0, 2.0, 3.0, 0.0, 3.0, 1.0];
833        let indices = vec![0, 0, 1, 2, 0, 1, 2];
834        let indptr = vec![0, 1, 3, 7];
835
836        let sym_matrix =
837            SymCsrMatrix::new(data, indptr, indices, (3, 3)).expect("Operation failed");
838        let sym_array = SymCsrArray::new(sym_matrix);
839
840        assert_eq!(sym_array.shape(), (3, 3));
841        assert!(is_sym_sparse_array(&sym_array));
842
843        // Check values
844        assert_eq!(SparseArray::get(&sym_array, 0, 0), 2.0);
845        assert_eq!(SparseArray::get(&sym_array, 0, 1), 1.0);
846        assert_eq!(SparseArray::get(&sym_array, 1, 0), 1.0); // Symmetric element
847        assert_eq!(SparseArray::get(&sym_array, 1, 2), 3.0);
848        assert_eq!(SparseArray::get(&sym_array, 2, 1), 3.0); // Symmetric element
849
850        // Convert to standard CSR
851        let csr = SymSparseArray::to_csr(&sym_array).expect("Operation failed");
852        assert_eq!(csr.nnz(), 10); // Full matrix with symmetric elements
853    }
854
855    #[test]
856    fn test_sym_coo_array() {
857        // Create a symmetric matrix in COO format
858        let data = vec![2.0, 1.0, 2.0, 3.0, 1.0];
859        let rows = vec![0, 1, 1, 2, 2];
860        let cols = vec![0, 0, 1, 1, 2];
861
862        let sym_matrix = SymCooMatrix::new(data, rows, cols, (3, 3)).expect("Operation failed");
863        let sym_array = SymCooArray::new(sym_matrix);
864
865        assert_eq!(sym_array.shape(), (3, 3));
866        assert!(is_sym_sparse_array(&sym_array));
867
868        // Check values
869        assert_eq!(SparseArray::get(&sym_array, 0, 0), 2.0);
870        assert_eq!(SparseArray::get(&sym_array, 0, 1), 1.0);
871        assert_eq!(SparseArray::get(&sym_array, 1, 0), 1.0); // Symmetric element
872        assert_eq!(SparseArray::get(&sym_array, 1, 2), 3.0);
873        assert_eq!(SparseArray::get(&sym_array, 2, 1), 3.0); // Symmetric element
874
875        // Test from_triplets with enforce symmetry
876        // Input is intentionally asymmetric - will be fixed by enforce_symmetric=true
877        let rows2 = vec![0, 0, 1, 1, 2, 1, 0];
878        let cols2 = vec![0, 1, 1, 2, 2, 0, 2];
879        let data2 = vec![2.0, 1.5, 2.0, 3.5, 1.0, 0.5, 0.0];
880
881        let sym_array2 = SymCooArray::from_triplets(&rows2, &cols2, &data2, (3, 3), true)
882            .expect("Operation failed");
883
884        // Should average the asymmetric values
885        assert_eq!(SparseArray::get(&sym_array2, 0, 1), 1.0); // Average of 1.5 and 0.5
886        assert_eq!(SparseArray::get(&sym_array2, 1, 0), 1.0); // Symmetric element
887        assert_eq!(SparseArray::get(&sym_array2, 0, 2), 0.0); // Zero element
888    }
889
890    #[test]
891    fn test_construct_sym_utils() {
892        // Test creating an identity matrix
893        let eye = construct_sym::eye_sym_array::<f64>(3, "csr").expect("Operation failed");
894
895        assert_eq!(eye.shape(), (3, 3));
896        assert_eq!(SparseArray::get(&*eye, 0, 0), 1.0);
897        assert_eq!(SparseArray::get(&*eye, 1, 1), 1.0);
898        assert_eq!(SparseArray::get(&*eye, 2, 2), 1.0);
899        assert_eq!(SparseArray::get(&*eye, 0, 1), 0.0);
900
901        // Test creating a tridiagonal matrix - with coo format since csr had issues
902        let diag = vec![2.0, 2.0, 2.0];
903        let offdiag = vec![1.0, 1.0];
904
905        let tri =
906            construct_sym::tridiagonal_sym_array(&diag, &offdiag, "coo").expect("Operation failed");
907
908        assert_eq!(tri.shape(), (3, 3));
909        assert_eq!(SparseArray::get(&*tri, 0, 0), 2.0); // Main diagonal
910        assert_eq!(SparseArray::get(&*tri, 1, 1), 2.0);
911        assert_eq!(SparseArray::get(&*tri, 2, 2), 2.0);
912        assert_eq!(SparseArray::get(&*tri, 0, 1), 1.0); // Off-diagonal
913        assert_eq!(SparseArray::get(&*tri, 1, 0), 1.0); // Symmetric element
914        assert_eq!(SparseArray::get(&*tri, 1, 2), 1.0);
915        assert_eq!(SparseArray::get(&*tri, 0, 2), 0.0); // Zero element
916
917        // Test creating a banded matrix
918        let diagonals = vec![
919            vec![2.0, 2.0, 2.0, 2.0, 2.0], // Main diagonal
920            vec![1.0, 1.0, 1.0, 1.0],      // First off-diagonal
921            vec![0.5, 0.5, 0.5],           // Second off-diagonal
922        ];
923
924        let band = construct_sym::banded_sym_array(&diagonals, 5, "csr").expect("Operation failed");
925
926        assert_eq!(band.shape(), (5, 5));
927        assert_eq!(SparseArray::get(&*band, 0, 0), 2.0);
928        assert_eq!(SparseArray::get(&*band, 0, 1), 1.0);
929        assert_eq!(SparseArray::get(&*band, 0, 2), 0.5);
930        assert_eq!(SparseArray::get(&*band, 2, 0), 0.5); // Symmetric element
931    }
932
933    #[test]
934    fn test_sym_conversions() {
935        // Create a symmetric matrix
936        // Lower triangular part only
937        let data = vec![2.0, 1.0, 2.0, 3.0, 1.0];
938        let rows = vec![0, 1, 1, 2, 2];
939        let cols = vec![0, 0, 1, 1, 2];
940
941        let sym_coo = SymCooArray::from_triplets(&rows, &cols, &data, (3, 3), true)
942            .expect("Operation failed");
943
944        // Convert to symmetric CSR
945        let sym_csr = sym_coo.to_sym_csr().expect("Operation failed");
946
947        // Check values are preserved
948        for i in 0..3 {
949            for j in 0..3 {
950                assert_eq!(
951                    SparseArray::get(&sym_coo, i, j),
952                    SparseArray::get(&sym_csr, i, j)
953                );
954            }
955        }
956
957        // Convert to standard formats
958        let csr = SymSparseArray::to_csr(&sym_coo).expect("Operation failed");
959        let coo = SymSparseArray::to_coo(&sym_csr).expect("Operation failed");
960
961        // Check full symmetric matrix in standard formats
962        assert_eq!(csr.nnz(), 7); // Accounts for symmetric pairs
963        assert_eq!(coo.nnz(), 7);
964
965        for i in 0..3 {
966            for j in 0..3 {
967                assert_eq!(SparseArray::get(&csr, i, j), SparseArray::get(&coo, i, j));
968                assert_eq!(
969                    SparseArray::get(&csr, i, j),
970                    SparseArray::get(&sym_csr, i, j)
971                );
972            }
973        }
974    }
975}
scirs2_sparse/lib.rs

scirs2_sparse/
lib.rs
