Allow Matrix diagonal as mutable elements about rulinalg HOT 19 CLOSED

I like this idea - though it is possible that it may be unsafe in places (due to assumptions in other areas of the code). It may be fine but we'd need to check carefully that none of our iterators etc. broke in any weird cases.

We should have assert! checks to ensure that the inputs are good, then use get_unchecked_mut if all is well.

I'm not too sure what you meant by this part though? If all worked as expected we should be able to access, iterate, etc. the elements of the slice in the same way as we would normally.

from rulinalg.

I was thinking about this too, the other day. I think we could perhaps take it one step further, by letting us take slices also from off-diagonals? I.e., we introduce some convention to number the diagonals, and then to get the kth diagonal, you would call something like mat.diagonal(k). Or maybe mat.diagonal() should be reserved for the standard diagonal, and mat.off_diagonal(k) would be better. I'm just thinking out loud here. MATLAB has the following convention, by the way:

D = diag(v,k) places the elements of vector v on the kth diagonal. k=0 represents the main diagonal, k>0 is above the main diagonal, and k<0 is below the main diagonal.

@AtheMathmo: In terms of checking carefully that none of our iterators break in weird cases, I think that's an excellent opportunity to write more tests for said iterators!

Edit: I have to say that I think MATLAB's convention is somewhat opposite of what I'd prefer, though.

from rulinalg.

Thinking of a function signature like:

impl Matrix ... 
    pub fn diag_slice_mut(&self mut, start_row: usize) -> MatrixSliceMut;

This would allow us to get diagonals aside from just the main diagonal. Basically what @Andlon just said. I'm not convinced we should follow MATLAB's convention, but let's look at a few other examples first and see what would be the most elegant for this library.

from rulinalg.

@andrewcsmith: Agree. Do note, however, that in your case it's impossible to slice above the main diagonal, since usize cannot be negative! Also, a negative start_row would be a bit strange, no...?

from rulinalg.

Right, I'm not tied to this either. Seems that (no surprise) numpy also uses k so maybe that's the way to go.

from rulinalg.

It doesn't play very nice with our current usize-based indexing though... Slightly worried we'll introduce weird corner cases from this. Any thoughts?

from rulinalg.

I do think this is really valuable functionality - and I think it is worth working on. But a few more thoughts/concerns we should address:

• We'd like to have this in BaseMatrix ideally. But how do we handle calling this/these functions on a diagonal slice?
• On a similar note to the above - I have pretty liberally assumed that a MatrixSlice is a contiguous (in the 2-d sense) block of data. This may break things too.
• We'd want to make clear that this diag has n rows and n cols (where n = min(rows, cols)). But having no. elements != n*n will probably cause some breakage.

If we did want to have a single function for getting the diagonal (which I think we do) - maybe an enum is an appropriate way to go?

enum Diag {
    Above(usize),
    Below(usize),
    [Main]?
}

(Forgive the naming)

from rulinalg.

I don't think that k is an index if we don't treat it like one. In the case of matlab etc. k is an offset, not an index.

I agree, slice is not the right word for this. In my opinion, slice should pretty much always mean contiguous. (At least the MatrixSlice has contiguous rows.) Instead, perhaps we should just use an iterator.

My main concern is making sure we have a way to iterate over the elements mutably, rather than iterate over the indices and use get_unchecked. If we can reduce the scope of the unsafe access and use an iterator instead it will be an improvement. Might be worthwhile adding some benchmarks to prevent regressions.

from rulinalg.

I think an Iterator is a better idea - mostly because it will be a lot easier to implement. Using the current MatrixSlice would be near impossible I think - all the arithmetic would also be prone to errors.

I think adding a Diagonal(Mut) iterator which can be used to iterate over the (off-)diagonal elements is a smart move. We can retain existing functionality by collecting the iterator - though we'll probably want to force a breaking change on the diag function (so it returns this iterator instead of a Vec). We can modify the from_diag function to take an Into<Diagonal> (and implement this for Vec/&[T]).

This seems a lot more attainable (and less frightening) to me.

from rulinalg.

@AtheMathmo: Definitely agree that an iterator is the way to go here!

from rulinalg.

By the way, one thing that's very useful about how Numpy and MATLAB models super/sub-diagonals, is that by having a single number, you can write code that, for example, scales elements in a three-diagonal matrix:

// Might not work exactly as I've written here, but I hope you get the point!
let mut diagonals = (-1 .. 1).flat_map(|k| mat.diag_mut(k));
for ref mut element in diagonals {
    *element *= scale_factor;
}

So there definitely is some merit to this convention.

from rulinalg.

I think I understand - though I'm not sure that the borrow checker would be happy with that? I think that multiple mutable references to mat would be an issue that we couldn't get around.

from rulinalg.

Ah, good catch! That's a shame. My example also has the wrong range (it will only include -1 and 0), by the way.

from rulinalg.

Anyway, I think my point still holds, as you can rewrite that particular example to make the borrow checker happy.

// The -1 .. 2 range doesn't look so intuitive though
for k in -1 .. 2 {
    for ref mut element in mat.diag_mut(k) {
        *element *= scale_factor;
    }
}

from rulinalg.

You can! It's possible that the compiler is smart enough to know that one borrow ends before the next begins in the flat_map too actually.

I'm pretty sold on the iterator - if no one else wants to I may try to get a PR in within the next 24 hours (or less). I'll stick with non-breaking changes for now though.

from rulinalg.

What, really? Wow, rustc never ceases to amaze! That's really cool, actually :) I'm all in favor!

from rulinalg.

Ah, sorry maybe I was unclear. I haven't tested it - I just thought it might be possible :) (and would indeed be amazing).

from rulinalg.

Oh, heh, all right. I thought perhaps you'd tried a simplified example with placeholders already. Well, in any case, I suppose my other proposed alternative is almost as good.

from rulinalg.

I have a WIP PR in place for this. There are still some questions to resolve that we discussed here.

from rulinalg.