Skip to content

Commit 99663dc

Browse files
committed
Prototype concrete indptr type
As suggested by @mulimoen (#244 (comment)), it's doable to have a concrete type abstracting the classic indptr case with the first element at 0, and the sliced matrix case, by simply performing the subtraction everywhere.
1 parent f3db918 commit 99663dc

File tree

2 files changed

+141
-51
lines changed

2 files changed

+141
-51
lines changed

src/lib.rs

Lines changed: 6 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -89,13 +89,12 @@ pub use crate::indexing::SpIndex;
8989

9090
pub use crate::sparse::{
9191
csmat::CsIter, csmat::OuterIterator, csmat::OuterIteratorMut,
92-
csmat::OuterIteratorPerm, indptr::IndPtrStorage,
93-
kronecker::kronecker_product, CsMat, CsMatBase, CsMatI, CsMatVecView,
94-
CsMatView, CsMatViewI, CsMatViewMut, CsMatViewMutI, CsStructure,
95-
CsStructureI, CsStructureView, CsStructureViewI, CsVec, CsVecBase, CsVecI,
96-
CsVecView, CsVecViewI, CsVecViewMut, CsVecViewMutI, SparseMat, TriMat,
97-
TriMatBase, TriMatI, TriMatIter, TriMatView, TriMatViewI, TriMatViewMut,
98-
TriMatViewMutI,
92+
csmat::OuterIteratorPerm, indptr::IndPtr, kronecker::kronecker_product,
93+
CsMat, CsMatBase, CsMatI, CsMatVecView, CsMatView, CsMatViewI,
94+
CsMatViewMut, CsMatViewMutI, CsStructure, CsStructureI, CsStructureView,
95+
CsStructureViewI, CsVec, CsVecBase, CsVecI, CsVecView, CsVecViewI,
96+
CsVecViewMut, CsVecViewMutI, SparseMat, TriMat, TriMatBase, TriMatI,
97+
TriMatIter, TriMatView, TriMatViewI, TriMatViewMut, TriMatViewMutI,
9998
};
10099

101100
pub use crate::sparse::symmetric::is_symmetric;

src/sparse/indptr.rs

Lines changed: 135 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -1,66 +1,157 @@
11
//! This module defines the behavior of types suitable to be used
22
//! as `indptr` storage in a [`CsMatBase`].
3-
///
4-
/// [`CsMatBase`]: type.CsMatBase.html
3+
//!
4+
//! [`CsMatBase`]: type.CsMatBase.html
5+
6+
use crate::errors::SprsError;
57
use crate::indexing::SpIndex;
8+
use std::ops::Deref;
69

7-
/// Iterator type for types implementing `IndPtrStorage`. Since the language
8-
/// does not allow `impl Iterator` in traits yet, we have to use a concrete
9-
/// type.
10-
pub enum IptrIter<'a, Iptr> {
11-
Trivial(std::slice::Windows<'a, Iptr>),
12-
Offset(Iptr, std::slice::Windows<'a, Iptr>),
10+
pub struct IndPtr<Iptr, Storage>
11+
where
12+
Iptr: SpIndex,
13+
Storage: Deref<Target = [Iptr]>,
14+
{
15+
storage: Storage,
1316
}
1417

15-
impl<'a, Iptr: SpIndex> Iterator for IptrIter<'a, Iptr> {
16-
type Item = (Iptr, Iptr);
17-
18-
fn next(&mut self) -> Option<<Self as Iterator>::Item> {
19-
match self {
20-
IptrIter::Trivial(iter) => iter.next().map(|x| (x[0], x[1])),
21-
IptrIter::Offset(offset, iter) => {
22-
iter.next().map(|x| (x[0] + *offset, x[1] + *offset))
18+
impl<Iptr, Storage> IndPtr<Iptr, Storage>
19+
where
20+
Iptr: SpIndex,
21+
Storage: Deref<Target = [Iptr]>,
22+
{
23+
pub(crate) fn check_structure(storage: &Storage) -> Result<(), SprsError> {
24+
for i in storage.iter() {
25+
if i.try_index().is_none() {
26+
return Err(SprsError::IllegalArguments(
27+
"Indptr value out of range of usize",
28+
));
2329
}
2430
}
31+
if !storage
32+
.windows(2)
33+
.all(|x| x[0].index_unchecked() <= x[1].index_unchecked())
34+
{
35+
return Err(SprsError::UnsortedIndptr);
36+
}
37+
if storage
38+
.last()
39+
.cloned()
40+
.map(Iptr::index_unchecked)
41+
.map(|i| i > usize::max_value() / 2)
42+
.unwrap_or(false)
43+
{
44+
// We do not allow indptr values to be larger than half
45+
// the maximum value of an usize, as that would clearly exhaust
46+
// all available memory
47+
// This means we could have an isize, but in practice it's
48+
// easier to work with usize for indexing.
49+
return Err(SprsError::IllegalArguments(
50+
"An indptr value is larger than allowed",
51+
));
52+
}
53+
Ok(())
54+
}
55+
56+
pub fn new(storage: Storage) -> Result<Self, crate::errors::SprsError> {
57+
IndPtr::check_structure(&storage).map(|_| IndPtr::new_trusted(storage))
2558
}
26-
}
2759

28-
/// Index Pointer Storage specification.
29-
///
30-
/// A type implementing this trait can be used in sparse matrix algorithms in
31-
/// the CSR or CSC format.
32-
pub trait IndPtrStorage<Iptr: SpIndex> {
33-
/// The length of the index pointer storage.
34-
fn len(&self) -> usize;
60+
pub(crate) fn new_trusted(storage: Storage) -> Self {
61+
IndPtr { storage }
62+
}
3563

36-
/// Access the index pointer element at location `i`. Returns `None`
37-
/// if `i >= self.len()`.
38-
fn get(&self, i: usize) -> Option<Iptr>;
64+
/// The length of the underlying storage
65+
pub fn len(&self) -> usize {
66+
self.storage.len()
67+
}
3968

40-
/// Access the index pointer element at location `i`.
41-
///
42-
/// # Panics
69+
/// The number of outer dimensions this indptr represents
70+
pub fn outer_dims(&self) -> usize {
71+
if self.storage.len() >= 1 {
72+
self.storage.len() - 1
73+
} else {
74+
0
75+
}
76+
}
77+
78+
/// Indicates whether the underlying storage is proper, which means the
79+
/// indptr corresponds to a non-sliced matrix.
4380
///
44-
/// If `i >= self.len()`
45-
fn at(&self, i: usize) -> Iptr {
46-
<Self as IndPtrStorage<Iptr>>::get(self, i).unwrap()
81+
/// An empty matrix is considered non-proper.
82+
pub fn is_proper(&self) -> bool {
83+
self.storage
84+
.get(0)
85+
.map(|i| *i == Iptr::zero())
86+
.unwrap_or(false)
4787
}
4888

49-
/// Iterate over elements in the index pointer storage
50-
fn iter(&self) -> IptrIter<Iptr>;
51-
}
89+
/// Return a view on the underlying slice if it is a proper `indptr` slice,
90+
/// which is the case if its first element is 0. `None` will be returned
91+
/// otherwise.
92+
pub fn slice(&self) -> Option<&[Iptr]> {
93+
if self.is_proper() {
94+
Some(&self.storage[..])
95+
} else {
96+
None
97+
}
98+
}
5299

53-
impl<'a, Iptr: SpIndex> IndPtrStorage<Iptr> for &'a [Iptr] {
54-
fn len(&self) -> usize {
55-
<[Iptr]>::len(self)
100+
/// Return a view of the underlying storage. Should be used with care in
101+
/// sparse algorithms, as this won't check if the storage corresponds to a
102+
/// proper matrix
103+
pub fn raw_storage(&self) -> &[Iptr] {
104+
&self.storage[..]
56105
}
57106

58-
fn get(&self, i: usize) -> Option<Iptr> {
59-
<[Iptr]>::get(self, i).map(|x| *x)
107+
fn offset(&self) -> Iptr {
108+
self.storage.get(0).cloned().unwrap_or(Iptr::zero())
109+
}
110+
111+
/// Iterate over outer dimensions, yielding start and end indices for each
112+
/// outer dimension.
113+
pub fn iter_outer(
114+
&self,
115+
) -> impl std::iter::DoubleEndedIterator<Item = (Iptr, Iptr)>
116+
+ std::iter::DoubleEndedIterator<Item = (Iptr, Iptr)>
117+
+ '_ {
118+
let offset = self.offset();
119+
self.storage.windows(2).map(move |x| {
120+
if offset == Iptr::zero() {
121+
(x[0], x[1])
122+
} else {
123+
(x[0] - offset, x[1] - offset)
124+
}
125+
})
126+
}
127+
128+
/// Return the value of the indptr at index i. This method is intended for
129+
/// low-level usage only, `outer_inds` should be preferred most of the time
130+
pub fn index(&self, i: usize) -> Iptr {
131+
let offset = self.offset();
132+
self.storage[i] - offset
133+
}
134+
135+
/// Get the start and end indices for the requested outer dimension
136+
///
137+
/// # Panics
138+
///
139+
/// If `i >= self.outer_dims()`
140+
pub fn outer_inds(&self, i: usize) -> (Iptr, Iptr) {
141+
assert!(i + 1 < self.storage.len());
142+
let offset = self.offset();
143+
(self.storage[i] - offset, self.storage[i + 1] - offset)
60144
}
61145

62-
/// Iterate over elements in the index pointer storage
63-
fn iter(&self) -> IptrIter<Iptr> {
64-
IptrIter::Trivial(self.windows(2))
146+
/// The number of nonzero elements described by this indptr
147+
pub fn nnz(&self) -> usize {
148+
let offset = self.offset();
149+
// index_unchecked validity: structure checks ensure that the last index
150+
// larger than the first, and that both can be represented as an usize
151+
self.storage
152+
.last()
153+
.map(|i| *i - offset)
154+
.map(Iptr::index_unchecked)
155+
.unwrap_or(0)
65156
}
66157
}

0 commit comments

Comments
 (0)