Matrix.hpp

Go to the documentation of this file.

#ifndef RAIIGRAPH_MATRIX_HPP
#define RAIIGRAPH_MATRIX_HPP
 
#include "igraph.h"
#include "error.hpp"
#include "Vector.hpp"
#include <algorithm>
#include <iterator>
 
namespace raiigraph {
 
template<class Ns_>
class Matrix {
private:
    void setup(igraph_integer_t nr, igraph_integer_t nc) {
        if (Ns_::init(&my_matrix, nr, nc)) {
            throw std::runtime_error("failed to initialize igraph matrix of dimensions " + std::to_string(nr) + " x " + std::to_string(nc));
        }
    }
 
public:
    typedef typename Ns_::igraph_type igraph_type;
 
    typedef typename Ns_::value_type value_type;
 
    typedef value_type& reference;
 
    typedef const value_type& const_reference;
 
    typedef igraph_integer_t size_type;
 
    typedef igraph_integer_t difference_type;
 
    typedef value_type* iterator;
 
    typedef const value_type* const_iterator;
 
    typedef std::reverse_iterator<iterator> reverse_iterator;
 
    typedef std::reverse_iterator<const_iterator> reverse_const_iterator;
 
    typedef typename Ns_::vector_type vector_type;
 
public:
    Matrix() : Matrix(0, 0) {}
 
    Matrix(size_type nr, size_type nc, const value_type& val = value_type()) {
        setup(nr, nc);
        if (val != 0) { // setup() already zero-initializes the backing array.
            std::fill(begin(), end(), val);
        }
    }
 
    Matrix(igraph_type&& matrix) : my_matrix(std::move(matrix)) {}
 
public:
    Matrix(const Matrix<Ns_>& other) {
        if (Ns_::copy(&my_matrix, &(other.my_matrix))) {
            throw std::runtime_error("failed to copy-construct integer igraph matrix");
        }
    }
 
    Matrix<Ns_>& operator=(const Matrix<Ns_>& other) {
        if (this != &other) {
            if (Ns_::update(&my_matrix, &(other.my_matrix))) { // my_matrix should already be initialized before the assignment.
                throw std::runtime_error("failed to copy-assign integer igraph matrix");
            }
        }
        return *this;
    }
 
    Matrix(Matrix<Ns_>&& other) {
        setup(0, 0); // we must leave 'other' in a valid state.
        std::swap(my_matrix, other.my_matrix);
    }
 
    Matrix& operator=(Matrix<Ns_>&& other) {
        if (this != &other) {
            std::swap(my_matrix, other.my_matrix); // 'my_matrix' should already be initialized, so we're leaving 'other' in a valid state.
        }
        return *this;
    }
 
    ~Matrix() {
        Ns_::destroy(&my_matrix);
    }
 
public:
    igraph_bool_t empty() const {
        return Ns_::empty(&my_matrix);
    }
 
    size_type size() const {
        return Ns_::size(&my_matrix);
    }
 
    size_type nrow() const {
        return my_matrix.nrow; 
    }
 
    size_type ncol() const {
        return my_matrix.ncol;
    }
 
    constexpr size_type max_size() const { return IGRAPH_INTEGER_MAX; }
 
    size_type capacity() const {
        return my_matrix.data.stor_end - my_matrix.data.stor_begin;
    }
 
public:
    void clear() {
        resize(0, 0);
    }
 
    void resize(size_type nr, size_type nc, value_type val = value_type()) {
        auto old_size = this->size();
        check_code(Ns_::resize(&my_matrix, nr, nc));
        auto new_size = this->size();
        if (old_size < new_size) {
            std::fill_n(begin() + old_size, new_size - old_size, val);
        }
    }
 
    void shrink_to_fit() {
        Ns_::shrink_to_fit(&my_matrix);
    }
 
public:
    reference operator[](size_type i) {
        return *(begin() + i);
    }
 
    const_reference operator[](size_type i) const {
        return *(begin() + i);
    }
 
    reference operator()(size_type r, size_type c) {
        return *(begin() + r + c * my_matrix.nrow); // no need to worry about overflow as size_type is guaranteed to hold the full size.
    }
 
    const_reference operator()(size_type r, size_type c) const {
        return *(begin() + r + c * my_matrix.nrow);
    }
 
    reference back() {
        return *(end() - 1);
    }
 
    const_reference back() const {
        return *(end() - 1);
    }
 
    reference front() {
        return *(begin());
    }
 
    const_reference front() const {
        return *(begin());
    }
 
public:
    iterator begin() {
        return my_matrix.data.stor_begin;
    }
 
    iterator end() {
        return my_matrix.data.stor_end;
    }
 
    const_iterator begin() const {
        return cbegin();
    }
 
    const_iterator end() const {
        return cend();
    }
 
    const_iterator cbegin() const {
        return my_matrix.data.stor_begin;
    }
 
    const_iterator cend() const {
        return my_matrix.data.stor_end;
    }
 
    value_type* data() {
        return my_matrix.data.stor_begin;
    }
 
    const value_type* data() const {
        return my_matrix.data.stor_begin;
    }
 
    reverse_iterator rbegin() {
        return std::reverse_iterator(end());
    }
 
    reverse_iterator rend() {
        return std::reverse_iterator(begin());
    }
 
    reverse_const_iterator rbegin() const {
        return std::reverse_iterator(end());
    }
 
    reverse_const_iterator rend() const {
        return std::reverse_iterator(begin());
    }
 
    reverse_const_iterator crbegin() const {
        return std::reverse_iterator(cend());
    }
 
    reverse_const_iterator crend() const {
        return std::reverse_iterator(cbegin());
    }
 
public:
    template<typename BaseIterator, typename BaseReference>
    class View {
    public:
        View(BaseIterator start, size_type step_size, size_type max_steps) : start(start), step_size(step_size), max_steps(max_steps) {}
 
    private:
        BaseIterator start;
        size_type step_size, max_steps;
    public:
        bool empty() const {
            return max_steps == 0;
        }
 
        size_type size() const {
            return max_steps;
        }
 
        BaseReference operator[](size_type i) const {
            return *(start + i * step_size); // no need to cast to avoid overflow, as size_type determines the max size anyway.
        }
 
        BaseReference back() const {
            return (*this)[max_steps - 1];
        }
 
        BaseReference front() const {
            return *start;
        }
 
    public:
        struct Iterator {
        private:
            BaseIterator start;
            size_type step_size = 0;
 
            // Note that we don't just shift 'start' directly as defining the
            // 'end()' iterator of the view could up shifting 'start' to an
            // address beyond the one-past-the-end of the backing array in the
            // 'igraph_matrix_t'. Pointer arithmetic beyond the array bounds is
            // undefined IIRC, so instead, we store the offset and att it to
            // 'start' upon dereference of the iterator. This ensures that an
            // invalid address is never constructed even if it is unused.
            size_type offset = 0; 
 
        public:
            using iterator_category = std::random_access_iterator_tag;
            typedef typename Ns_::value_type value_type;
            typedef decltype(&std::declval<BaseReference>()) pointer;
            typedef igraph_integer_t difference_type;
            typedef BaseReference reference;
 
        public:
            explicit Iterator(BaseIterator start, size_type step_size, size_type position) : 
                start(std::move(start)), step_size(step_size), offset(position * step_size) {}
 
            Iterator() = default;
 
        public:
            // List of required methods taken from https://cplusplus.com/reference/iterator/RandomAccessIterator/.
            bool operator==(const Iterator& other) const { 
                // see comments at https://stackoverflow.com/questions/4657513/comparing-iterators-from-different-containers
                // regarding the UB of comparing iterators from different containers, so we'll just quietly ignore it.
                return offset == other.offset;  
            }
 
            bool operator!=(const Iterator& other) const { 
                return offset != other.offset;
            }
 
            bool operator<(const Iterator& other) const {
                return offset < other.offset;
            }
 
            bool operator>(const Iterator& other) const {
                return offset > other.offset;
            }
 
            bool operator<=(const Iterator& other) const {
                return offset <= other.offset;
            }
 
            bool operator>=(const Iterator& other) const {
                return offset >= other.offset;
            }
 
        public:
            BaseReference operator*() const {
                return *(start + offset);
            }
 
            BaseReference operator[](size_type i) const {
                return *(start + offset + step_size * i);
            }
 
        public:
            Iterator& operator++() { 
                offset += step_size;
                return *this; 
            }
 
            Iterator operator++(int) { 
                auto copy = *this;
                offset += step_size;
                return copy;
            }
 
            Iterator& operator--() { 
                offset -= step_size;
                return *this; 
            }
 
            Iterator operator--(int) { 
                auto copy = *this;
                offset -= step_size;
                return copy;
            }
 
            Iterator operator+(difference_type n) const { 
                auto copy = *this;
                copy.offset += step_size * n;
                return copy;
            }
 
            Iterator& operator+=(difference_type n) { 
                offset += step_size * n;
                return *this; 
            }
 
            friend Iterator operator+(difference_type n, Iterator it) {
                it.offset += it.step_size * n;
                return it;
            }
 
            Iterator operator-(difference_type n) const { 
                auto copy = *this;
                copy.offset -= step_size * n;
                return copy;
            }
 
            Iterator& operator-=(difference_type n) { 
                offset -= step_size * n;
                return *this;
            }
 
            difference_type operator-(const Iterator& other) const {
                if (other.offset > offset) {
                    difference_type delta = (other.offset - offset) / step_size;
                    return -delta;
                } else {
                    return (offset - other.offset) / step_size;
                }
            }
        };
 
        Iterator begin() const {
            return Iterator(start, step_size, 0);
        }
 
        Iterator end() const {
            return Iterator(start, step_size, max_steps);
        }
 
        Iterator cbegin() const {
            return begin();
        }
 
        Iterator cend() const {
            return end();
        }
 
        std::reverse_iterator<Iterator> rbegin() const {
            return std::reverse_iterator(end());
        }
 
        std::reverse_iterator<Iterator> rend() const {
            return std::reverse_iterator(begin());
        }
 
        std::reverse_iterator<Iterator> crbegin() const {
            return std::reverse_iterator(cend());
        }
 
        std::reverse_iterator<Iterator> crend() const {
            return std::reverse_iterator(cbegin());
        }
    };
 
    View<iterator, reference> row(size_type r) {
        return View<iterator, reference>(begin() + r, my_matrix.nrow, my_matrix.ncol);
    }
 
    View<const_iterator, const_reference> row(size_type r) const {
        return View<const_iterator, const_reference>(begin() + r, my_matrix.nrow, my_matrix.ncol);
    }
 
    vector_type row_copy(size_type r) const {
        vector_type output(ncol());
        check_code(Ns_::get_row(&my_matrix, output.get(), r));
        return output;
    }
 
    View<iterator, reference> column(size_type c) {
        return View<iterator, reference>(begin() + c * my_matrix.nrow, 1, my_matrix.nrow);
    }
 
    View<const_iterator, const_reference> column(size_type c) const {
        return View<const_iterator, const_reference>(begin() + c * my_matrix.nrow, 1, my_matrix.nrow);
    }
 
    vector_type column_copy(size_type c) const {
        vector_type output(nrow());
        check_code(Ns_::get_col(&my_matrix, output.get(), c));
        return output;
    }
 
public:
    operator igraph_type*() {
        return &my_matrix;
    }
 
    operator const igraph_type*() const {
        return &my_matrix;
    }
 
    igraph_type* get() {
        return &my_matrix;
    }
 
    const igraph_type* get() const {
        return &my_matrix;
    }
 
public:
    void swap(Matrix<Ns_>& other)  {
        // Swapping structures entirely to ensure that iterators and pointers
        // remain valid; looks like igraph_matrix_swap does the same.
        std::swap(my_matrix, other.my_matrix);
    }
 
private:
    igraph_type my_matrix;
};
 
namespace matrix_internal {
 
struct Integer {
    typedef igraph_integer_t value_type;
    typedef igraph_matrix_int_t igraph_type;
    typedef IntVector vector_type;
 
#define RAIIGRAPH_MATRIX_SUFFIX _int
#include "fragments/matrix.hpp"
#undef RAIIGRAPH_MATRIX_SUFFIX
};
 
struct Real {
    typedef igraph_real_t value_type;
    typedef igraph_matrix_t igraph_type;
    typedef RealVector vector_type;
 
#define RAIIGRAPH_MATRIX_SUFFIX
#include "fragments/matrix.hpp"
#undef RAIIGRAPH_MATRIX_SUFFIX
};
 
struct Bool {
    typedef igraph_bool_t value_type;
    typedef igraph_matrix_bool_t igraph_type;
    typedef BoolVector vector_type;
 
#define RAIIGRAPH_MATRIX_SUFFIX _bool
#include "fragments/matrix.hpp"
#undef RAIIGRAPH_MATRIX_SUFFIX
};
 
}
typedef Matrix<matrix_internal::Integer> IntMatrix;
 
// For back-compatibility.
typedef IntMatrix IntegerMatrix;
typedef Matrix<matrix_internal::Real> RealMatrix;
 
typedef Matrix<matrix_internal::Bool> BoolMatrix;
 
}
 
 
#endif