sparsemat / sprs

sparse linear algebra library for rust
Apache License 2.0
381 stars 45 forks source link

Draft: Support generic dense vectors in more APIs #272

Closed vbarrielle closed 3 years ago

vbarrielle commented 3 years ago

Lots of APIs were using slices as their input, particularly when using output buffers. However, this was not a good fit for linear algebra, since consumers of the APIs are more likely to be using eg an array from ndarray.

Here we extend the DenseVector trait with a mutable version to be able to use it on output parameters. Since these are sealed traits we should be able to add unsafe indexing if necessary without a breaking change. It should also be possible to support eg nalgebra in the future without too much trouble.

This should improve the situation discussed in #93, though it's probably not done yet.

vbarrielle commented 3 years ago

This will probably require some iteration. @mulimoen I'd be very interested in your feedback :)

mulimoen commented 3 years ago

This looks great! I tried implementing DenseVector on Deref<Target = [N]> to replace some of the boilerplate implementations, but this seems to clash with ndarray::ArrayBase, even though this does not implement Deref. It seems rust has some forward protection based on the error message:

note: upstream crates may add a new impl of trait `std::ops::Deref` for type `ndarray::ArrayBase<_, ndarray::Dim<[usize; 1]>>` in future versions
mulimoen commented 3 years ago

I was thinking about incorporating the nshare crate, but I don't think this is a good idea. Adding support for nalgebra in the new traits is easy, and we would know there is no overhead. Outside crates should probably implement DenseVector themselves, without going through nshare.

vbarrielle commented 3 years ago

This looks great! I tried implementing DenseVector on Deref<Target = [N]> to replace some of the boilerplate implementations, but this seems to clash with ndarray::ArrayBase, even though this does not implement Deref.

Yes orphan rules will prevent us from having too broad implementations here, I think we need to implement on a case by case. Fortunately most implementations are trivial.

I was thinking about incorporating the nshare crate, but I don't think this is a good idea. Adding support for nalgebra in the new traits is easy, and we would know there is no overhead. Outside crates should probably implement DenseVector themselves, without going through nshare.

Right now the DenseVector trait is not publicly exposed by sprs, it is sealed. Which is a nice way for us to control precisely the types we want to support. I do think we want to support nalgebra better in the future, though probably with a feature flag to control the support.

I haven't played with nshare yet, so I don't know if it's convenient to use.

vbarrielle commented 3 years ago

With the latest changes, I think I'm leaning towards an interesting implementation, I've been able to write a single product impl allowing to apply a permutation to any dense vector (so a Vec, a slice, or an ndarray::Array).

I'm starting to think the DenseVector trait should be exposed outside the crate, I feel the need for it to make the solve in sprs-ldl more generic:

    //I'd like to use DenseVector here, but I probably need to expose it for public use
    pub fn solve<'a, V>(&self, rhs: &V) -> Vec<N>
    where
        N: 'a + Copy + Num,
        V: Deref<Target = [N]>,
    {
        let mut x = &self.symbolic.perm * &rhs[..];
        let l = self.l();
        ldl_lsolve(&l, &mut x);
        linalg::diag_solve(&self.diag, &mut x);
        ldl_ltsolve(&l, &mut x);
        let pinv = self.symbolic.perm.inv();
        &pinv * &x
    }

There's another interesting opportunity with this refactor: if I manage to implement everything that's related to dense vectors using these traits, then I can make the ndarray dependency optional (but on by default I guess), and also add an optional dependency to nalgebra that would also implement this trait. Which is probably not enough, as I'd need also a trait to represent dense matrices to convert from sparse to dense, but it looks like an interesting future direction.

vbarrielle commented 3 years ago

So I think I've reached a good state, I've been able to use DenseVector in all required places. @mulimoen if you're OK with the current state then I'll merge.

vbarrielle commented 3 years ago

Thanks for the review @mulimoen !