Specialization for BlasMatrix or BlasSubMatrix in solutions

[ClementPernet]

Presently, in every solution, there is a specialization for BlasMatrix but none for BlasSubMatrix As a consequence, when calling e.g. det on submatrix, the fallback used is the templated Matrix code, with BlasElimination method, which allocates and builds a fresh BlasMatrix out of our BlasSubMatrix.

We could avoid this non-sense, by enforcing that every solution has a specialization for BlasSubMatrix and none for BlasMatrix. Indeed a BlasMatrix can be automatically casted into a BlasSubMatrix since there are constructors for this purpose.

[bdsaunders]

I like this. It sounds right to me. I have one doubt vis a vis automatic casting: It would meet my desire except that we have a custom of passing by reference and following that custom we couldn't get the benefit of an automatically created temporary.

On Mon, Sep 17, 2018 at 12:10 PM, Clément Pernet [email protected] wrote:

Presently, in every solution, there is a specialization for BlasMatrix but none for BlasSubMatrix As a consequence, when calling e.g. det on submatrix, the fallback used is the templated Matrix code, with BlasElimination method, which allocates and builds a fresh BlasMatrix out of our BlasSubMatrix.

We could avoid this non-sense, by enforcing that every solution has a specialization for BlasSubMatrix and none for BlasMatrix. Indeed a BlasMatrix can be automatically casted into a BlasSubMatrix since there are constructors for this purpose.

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/linbox-team/linbox/issues/136, or mute the thread https://github.com/notifications/unsubscribe-auth/ADk6I7S2ljkHix6AfkXli3dm2jCqXvckks5ub8lcgaJpZM4WsPW2 .

[Breush]

@bdsaunders As you proposed in https://github.com/linbox-team/linbox/pull/134#issuecomment-422027791, BlasMatrix could inherit from BlasSubMatrix, I don't see any issues with that and that would prevent that temporary allocation.

[pascalgiorgi]

I would not be so optimistic. I see at least a potential problem in such a situation. Indeed, in the new design of BlasSubmatrix we can now allow the use of constness. This means that when we make a submatrix of a const matrix the submatrix cannot modify its data. This was not the case before, in fact it was not handled properly.

To achieve such a feature I need to transfer the constness of the given matrix associated to its type into the templated type of the submatrix construction, i.e.

const BlasMatrix<X,Y> A; BlasSubmatrix<const BlasMatrix<X,Y> > B(A,0,0,k,k);

using BlasSubmatrix<BlasMatrix<X,Y> > B(A,0,0,k,k) is incorrect now.

There, if one want to make BlasMatrix to inherit from BlasSubmatrix I think we lose the constness feature.

[bdsaunders]

On reflection, I'd say the custom of passing matrices by reference is not a strong enough reason against the casting approach Clement suggests. The user experience will be identical and a brief comment will explain it to the developer to avoid excessive head scratching. However, I need to be convinced by implementation that the compiler can actually find among the templates the right call...

On the other hand, there may be an elegant enough solution to Pascal's const-ness concern about inheritance... This deserves more thought.

By way of comparison, we have a partial implementation of a pure virtual blackbox base class in bb.h, intended to serve as the only type needed to call solutions using blackbox methods. It has it's issues including that read() and particularly rebind are problematic. I have sketched solutions in my head, but they are not very pretty. (template member classes and functions do not play well with pure virtuality.). It is possible this could be replaced by a light BB type to which other matrix reps are cast. It seems the issues about rebind remain...

On Tue, Sep 18, 2018 at 3:49 AM, Pascal Giorgi [email protected] wrote:

I would not be so optimistic. I see at least a potential problem in such a situation. Indeed, in the new design of BlasSubmatrix we can now allow the use of constness. This means that when we make a submatrix of a const matrix the submatrix cannot modify its data. This was not the case before, in fact it was not handled properly.

To achieve such a feature I need to transfer the constness of the given matrix associated to its type into the templated type of the submatrix construction, i.e.

const BlasMatrix<X,Y> A; BlasSubmatrix<const BlasMatrix<X,Y> > B(A,0,0,k,k);

using BlasSubmatrix<BlasMatrix<X,Y> > B(A,0,0,k,k) is incorrect now.

There, if one want to make BlasMatrix to inherit from BlasSubmatrix I think we lose the constness feature.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/linbox-team/linbox/issues/136#issuecomment-422291025, or mute the thread https://github.com/notifications/unsubscribe-auth/ADk6Iwo8wHu0ghx_MY9-tcIYZQazJdiOks5ucKVtgaJpZM4WsPW2 .

[Breush]

I was just thinking of using shared_ptr having just one class.

Demo:

#include <iostream>
#include <memory>
#include <vector>

class Matrix {
public:
    Matrix() {}

    // Basic constructor, creates the data
    Matrix(uint32_t m, uint32_t n)
        : _m(m)
        , _n(n)
        , _stride(n)
    {
        _data = std::make_shared<std::vector<float>>(m * n);
    }

    // Constructor making a sub matrix
    Matrix(const Matrix& M, uint32_t rowOffset, uint32_t columnOffset, uint32_t m, uint32_t n)
        : _data(M._data)
        , _m(m)
        , _n(n)
        , _startOffset(M._startOffset + rowOffset * M._stride + columnOffset)
        , _stride(M._stride)
    {
    }

    // operator= would share the data
    // copy() would really duplicate the data

    const float* getPointer() const { return _data->data() + _startOffset; }
    float* getWritePointer() { return _data->data() + _startOffset; }
    uint32_t stride() const { return _stride; }
    uint32_t rowdim() const { return _m; }
    uint32_t coldim() const { return _n; }

private:
    std::shared_ptr<std::vector<float>> _data;
    uint32_t _m = 0;
    uint32_t _n = 0;
    uint32_t _startOffset = 0;
    uint32_t _stride = 0;
};

#define MATRIX_DIM 16

int main(void)
{
    // We will store a matrix
    Matrix A;

    {
        // Creating data in an another scope
        Matrix M(MATRIX_DIM, MATRIX_DIM);
        for (auto i = 0u; i < MATRIX_DIM * MATRIX_DIM; ++i) {
            M.getWritePointer()[i] = i;
        }

        // Creating a submatrix
        A = Matrix(M, 1, 1, MATRIX_DIM - 2, MATRIX_DIM - 2);

        // M._data survives thanks to A._data
    }

    // A is writable!
    A.getWritePointer()[10] = 100;

    // Logging to be sure
    for (auto i = 0u; i < A.rowdim(); ++i) {
        for (auto j = 0u; j < A.coldim(); ++j) {
            std::cout << A.getPointer()[i * A.stride() + j] << " ";
        }
        std::cout << std::endl;
    }

    // A._data is destroyed, it's freed
    return 0;
}

It's a reference counting pointer, but I do not see what could go wrong. Any ideas?

[pascalgiorgi]

I am not sure that mixing matrix and submatrix in one class is a good idea. Using smart pointer could work but I don't want to pay the traversal of two pointers to access an element of the matrix.

Furthermore, we still have the problem of const-ness. Typically your example is misleading since your making a submatrix from a const matrix. Your program is working since you do not modify the data of the submatrix but it will not compile if you try to do A.setEntry(0,0,...). Basically, a submatrix must hold different type according to the const-ness of the original matrix.

[Breush]

I updated my example. Fact is the matrix data is writable even if we constructed from a const one because it is only one class and there is no restrictions. (But I guess friend classes would do the same?)

As you noticed, the bad thing is that shared_ptr to arrays directly are only in C++17, so I used a std::vector which indeed makes it dereference twice. However, I guess we can make our own "shared_ptr to a dynamic array" class.

[pascalgiorgi]

In fact that's my point. If I construct a submatrix from a const matrix I don't want that any data (even through pointer) are possibly modified, even if C++ allows this.

[Breush]

I see the point, I misunderstood. Still wanted to share my idea :)

[Breush]

Updated proposal for dense matrices:

• https://gist.github.com/Breush/a9802e59f1039a797436111cb2d11111

With in mind:

• Matrix<const Field> does not let the user modify data
• Matrix<Field> does let him do so
• A read-only submatrix of A is simply Matrix<const Field> B(A, rowOffset, colOffset, rowdim, coldim)
• You can do writable submatrices only if it comes from a writable one
• Submatrices keep data alive even if the parent is destroyed (thanks to shared_ptr)
• Only one pointer derefencing, as before
• I implemented operator[] just for fun

Pros:

• Only one class, no need for a submatrix specific one, and all computing routines (like det) can just implement Matrix<const Field> (as the non const version silently casts without trouble)
• The matrix data is kept alive as long as needed

Cons:

• <const Field> can feel a bit strange,

@pascalgiorgi Does that meet all criteria?

[ClementPernet]

I like the idea. Not sure I completely see how to adapt it with the design of container and iterators, but let see. Please change your uint32_t to size_t ;)