# Custom matrix class

It took me while to finish, but here is my custom matrix class. I assume that the row/column iterators are the most critical part of this class but anyway I would very much appreciate your ideas of this project.

#include <iostream>
#include <iomanip>
#include <ios>
#include <sstream>
#include <algorithm>
#include <functional>
#include <boost/range/algorithm.hpp>
#include <boost/typeof/typeof.hpp>
#include "Vector.h"

#ifndef MATRIX_H
#define MATRIX_H

template<class T,
std::size_t rowsize,
std::size_t colsize>
class Matrix{

private:

// row and/or colsize must be greater then 0;
static_assert(rowsize > 0, "Number of Rows must be greater then 0.");
static_assert(colsize > 0, "Number of Columns must be greater then 0.");

// Data is stored in a MVector, a modified std::vector
MVector<T> matrix;

// size of row dimension of the matrix
std::size_t row_dim;

// size of row dimension of the matrix
std::size_t column_dim;

public:

// *******************************************************
// Constructors
// *******************************************************

// default constructor
// @ param "rowsize" - row dimension of matrix
// @ param "colsize" - column dimension of matrix
Matrix() :matrix(rowsize * colsize),
row_dim(rowsize),
column_dim(colsize) {}

// initializer list constructor
// allows to create a matrix by an initialiser list
// example: Matrix<int, 2, 2> a_matrix = {1,2,3,4};
// @param "il" - initializer list
Matrix(std::initializer_list<T> il) :matrix(il),
row_dim(rowsize),
column_dim(colsize){}

// constructor by MVector
// if the size of the vector is larger then the size of the matrix,
// the matrix will construct with the first rowsize*colsize elements of the vector
// if the size of the vector is smaller then the size of the matrix,
// the matrix will not construct
explicit Matrix(const MVector<T>& s): matrix(s),
row_dim(rowsize),
column_dim(colsize){
}

explicit Matrix(MVector<T>&& s): matrix(s),
row_dim(rowsize),
column_dim(colsize){
}

// parameter pack constructor
// constructs a matrix by typing values into the brackets
// example: Matrix<int, 2, 2> a_matrix(1,2,3,4)
// values can be of any type, they will be interpreted as of class T
template<class ... N>
explicit Matrix(T first, N&&... values): matrix{first,
std::forward<T>(static_cast<T>(values))...},
row_dim(rowsize),
column_dim(colsize){}

// *******************************************************
// Iterator support
// *******************************************************

// traverse the entire vector
typename std::vector<T>::iterator Begin(){
return matrix.Begin();
}

typename std::vector<T>::iterator End(){
return matrix.End();
}

typename std::vector<T>::const_iterator Cbegin() const{
return matrix.Cbegin();
}

typename std::vector<T>::const_iterator Cend() const{
return matrix.Cend();
}

typename std::vector<T>::reverse_iterator Rbegin(){
return matrix.Rbegin();
}

typename std::vector<T>::reverse_iterator Rend(){
return matrix.Rend();
}

typename std::vector<T>::const_reverse_iterator Crbegin() const{
return matrix.Crbegin();
}

typename std::vector<T>::const_reverse_iterator Crend() const{
return matrix.Crend();
}

// column traversing
// don't unfortunatly lnow what to replace auto with
auto Begin_col( std::size_t i ){
}

auto End_col( std::size_t i ){
}

//row traversing
auto Begin_row( std::size_t i ){
}

auto End_row( std::size_t i ){
}

// returns the number of elements in the matrix
std::size_t size() const{
return matrix.size();
}

// returns the number of rows of the matrix
std::size_t rows() const{
return row_dim;
}

// returns the number of colums of the matrix
std::size_t cols() const{
return column_dim;
}

// returns the matrix as a vector
std::vector<T>& as_vector() {
return matrix.get_data();
}

// allows to access an element of the matrix by index expressed
// in terms of rows and columns
// @ param "r" - r'th row of the matrix
// @ param "c" - c'th column of the matrix
std::size_t index(std::size_t r, std::size_t c) const {
return r*cols()+c;
}

// returns a MVector object with the r'th row as it's elements
// slicing is possible from both ends and by "jumping" over elements
// @ param "begin" - starts at the n'th element
// @ param "end" - substracts m from from the last element.
// @ param "by" - selects every n'th element
MVector<T> get_row(std::size_t r, std::size_t begin = 0,
std::size_t end = 0, std::size_t by = 1) const{
MVector<T> row;
for (std::size_t i = index(r,begin); i < index(r,cols()-end); i += by) {
}
return row;
}

// get c'th column
// slicing is possible from both ends and by "jumping" over elements
// @ param "begin" - starts at the n'th element
// @ param "end" - substracts m from from the last element.
// @ param "by" - selects every n'th colum
MVector<T> get_column(std::size_t c, std::size_t begin = 0,
std::size_t end = 0, std::size_t by = 1) const{
assert(c < cols() && end < rows());
MVector<T> columns;
for (std::size_t i = index(begin, c); i < index(rows()-end,c); i+=by*cols()) {
}
return columns;
}

// get the diagonal elements of the matrix
// stors the diagonal in MVector
// and returns the vector
MVector<T> get_diagonal(){
MVector<T> diag;
Matrix<T, rowsize, colsize> temp = *this;
for (std::size_t i = 0; i < rowsize; i++) {
}
return diag;
}

// assignment operator
// assignes the content of another Matrix object to this one
// takes advantage of the overloaded " = " operator of the MVector class
// @ param "rhs" = Matrix object on the right hand side of the " = " sign
Matrix<T, rowsize, colsize>& operator=(Matrix<T, rowsize, colsize> rhs) {
this->matrix = rhs.matrix;
return *this;
}

// brackets operator
// return an elements stored in the matrix
// @ param "i" - i'th element in the matrix
T& operator[](std::size_t i) {
assert(i < size());
return matrix[i];
}

const T& operator[](std::size_t i) const {
assert(i < size());
return matrix[i];
}

// brackets operator
// return an elements stored in the matrix
// @ param "r" - r'th row in the matrix
// @ param "c" - c'th column in the matrix
T& operator()(std::size_t r, std::size_t c) {
assert(r < rows() && c < matrix.size() / rows());
return matrix[index(r,c)];
}

const T& operator()(std::size_t r, std::size_t c) const {
assert(r < rows() && c < matrix.size() / rows());
return matrix[index(r,c)];
}

// end of class
};

// *****************************************************************************
// // Equality operators
// *****************************************************************************

template<class T, std::size_t rowsize, std::size_t colsize>
bool operator==(const Matrix<T, rowsize, colsize> &lhs,
const Matrix<T, rowsize, colsize> &rhs) {

// defined in MVector:
return lhs.as_vector() == rhs.as_vector();
}
template<class T, std::size_t rowsize, std::size_t colsize>
bool operator!=(const Matrix<T, rowsize, colsize> &lhs,
const Matrix<T, rowsize, colsize> &rhs) {

// defined in MVector:
return (! (lhs == rhs) );
}

// *****************************************************************************
// // OS/iS  operators
// *****************************************************************************

template<class T, std::size_t row_dim, std::size_t col_dim>
std::ostream& operator<<(std::ostream &os, const Matrix<T, row_dim, col_dim>& m){
double max = *std::max_element( m.Cbegin(), m.Cend() );
std::stringstream convert;
convert << max;
for(std::size_t i = 0; i < m.rows(); i++){
os << "[" << i << ",]" << " ";
for(std::size_t j = 0; j <m.cols(); j++){
os << std::right << std::setw( convert.str().size() + 1 ) << m(i,j) << " ";
} os << std::endl;
}
return os;
}

// stores the data rowise
template<class T, std::size_t rowsize, std::size_t colsize>
std::istream& operator>>(std::istream &is, Matrix<T, rowsize, colsize>& m){

std::cout << "Enter " << rowsize*colsize << " numbers to construct the matrix: " <<std::endl;
for(std::size_t i = 0; i < rowsize*colsize; i++){

is >> m[i];
}
return is;
}

//more +=, +, *=, * ... overloaded operators


The first rule of matrix classes in C++ is that you don't write them by yourself; there are already dozens around, often highly optimized with [smart] expression templates and other cryptic stuff: Boost.uBLAS, Eigen, Blaze, Gmm++, Armadillo, Blitz++, etc...

That said, I would lie to you if I told you that I never tried to reimplement my own matrix class. I think that I tried to create one at leat two or three times :)

## Static versus dynamic storage

First of all, what I find odd is that your Matrix class, whose size is known at compile time, is built on the top of a MVector that is itself built on the top of an std::vector. That means that the memory is allocated of the heap and that your class will probably be much slower than intended. Since the size is known at compile time, you could have buit your class on the top of an std::array instead, which has way less overhead and should allocate its memory on the stack.

## Useless variables

The variable row_dim and column_dim are useless. They are known at compile time and the information is already in Matrix's template parameters rowsize and colsize (it would be good to capitalize template parameters by the way). You could get rid of these variables and reimplement the functions using them this way:

// returns the number of rows of the matrix
constexpr std::size_t rows() const{
return rowsize;
}

// returns the number of colums of the matrix
constexpr std::size_t cols() const{
return colsize;
}


By the way, I would expect functions named cols and rows to return columns and rows, not their size; you way want to change their names. You could also implement size so that it is also available at compile time:

constexpr std::size_t size() const{
return cols() * rows();
}


## Hide the implementation

Your class currently shows too many implementation details in its inteface. For example:

typename std::vector<T>::iterator Begin(){
return matrix.Begin();
}


The fact that your class uses an std::vector is an implementation detail. You should have written type aliases in Matrix and in MVector to provide a better abstraction:

// In Matrix
using iterator = typename MVector<T>::iterator;
using const_iterator = typename MVector<T>::const_iterator;


That way, if you change std::vector to something else in MVector, you won't have to update Matrix accordingly. It would also be a good idea to add other common subtypes such as value_type, pointer, reference... as in the standard library containers. A lot of template code out there relies on such subtypes.

## Standard library style

As I just said, writing containers that look like those used in the standard library generally makes them compatible with more existing code. Therefore, your functions should be named as in the standard library: instead of Begin and End, uncapitalized begin and end would have been better alternatives (that's actually why I try to write my code with the same case as the one used in the standard library). The full-lowercase member begin and end will work with both C++11's std::begin and std::end and with the range-based for loop.

Also, you miss a const overload for Begin and End.

## This should be a matrix by default, not a vector

Currently, many of the main operations clearly consider the class as a vector and not as a matrix. Those main operations are (c)begin, (c)end, operator[] and somehow the constructors. For a "matrix", I would expect operator[] to return a row or a column (it depends whether the matrix is row-major or column-major, documenting that would help) and (c)begin/(c)end to help iterating over rows or columns.

Being able to treat the matrix as a vector is often useful, but it should be explicitly required by the user, not the default option. If you want an example of implementation, the last matrix class I wrote provided the methods fbegin and fend that returned flat_iterators. I also wrote a container view, flat, so that I could write this:

matrix<int, 5, 4> mat = { /* ... */ };
for (int& val: flat(mat))
{
val *= 5;
}

• Thank you Morwenn for your comments. Indeed I was looking into boost matrix and eigen and still am very tempted to use them. However I am a C++ beginner and I must say say this was an excellent exercise. Aside of that, I feel more comfortable to customize my own class when I feel like changing something (like the ostream). – Vincent Oct 31 '14 at 11:14

There's not a lot left after the characteristically excellent answer from Morwenn, but I noticed a few things that may also help.

First, I see a lot right with this code. You didn't use using namespace std;, your functions are actually documented with comments and the code generally seems well though-out and complete. Good for you!

Here are some things that could be improved:

## Provide a matching #endif

This is probably just a cut-and-paste error, but there is no matching #endif for the #ifndef MATRIX_H at the top of the file.

## Fix the return type for some functions

As the comment indicates, you already knew the return type wasn't quite right for a few of your methods:

// don't unfortunatly lnow what to replace auto with
auto Begin_col( std::size_t i ){


I believe that the answer is that it's a boost::range_iterator<T>::type

typename boost::range_iterator<T>::type Begin_col( std::size_t i ){


However, with that said, it might be worth considering whether this is really the type you want to use for this.

## Fix operator<<

I supplied the missing Vector.h file and then used this test code:

int main()
{
Matrix<std::string, 2, 3> m{
"Larry", "Moe", "Curly",
"Belgium", "Netherlands", "Luxembourg"
};
std::cout << m << '\n';
}


This fails because of the code in the operator<< method:

std::ostream& operator<<(std::ostream &os, const Matrix<T, row_dim, col_dim>& m){
// uh oh... this assumes that the elements are numeric and
// can be converted to double
double max = *std::max_element( m.Cbegin(), m.Cend() );


I see what you're trying to do there, setting up the width so that all of the columns are aligned, but the method you're using isn't quite right. One way to do this would be to find the maximum width of all elements as printed and then use that to set the column width.

   const size_t colwidth = std::accumulate(m.Cbegin(), m.Cend(),
0, getPrintedLength<T>);


The helper function is then defined like this:

template<class T>
size_t getPrintedLength(size_t currmax, const T &thing)
{
std::ostringstream in;
in << thing;
return std::max(in.str().size(), currmax);
}


Now the column sizing works for any type for which std::ostream& operator<< is defined.

## Separate prompt from input in istream operator>>

Right now, the istream operator>> function prints a prompt and asks for $n$ numbers, where $n$ is the number of elements in the Matrix. However, as shown above, they aren't necessarily numbers. Or, if they are numbers, they might be complex numbers.

## Consider changing the Begin_col and End_col iterators

What happens if I iterate from m.Begin_col(0) to m.End_col(1)? I suspect it's not going to work as you expect, but I could be wrong since the functionality depends on the missing Vector.h code.

• Thank you for your comments Mr. Edwards. I don't think that I will be able to implement your suggestions as well as those from Morwenn immediatly. There is one thing about the std::array though. I have around 2 Gb of free DRAM on my computer. Thus I am able to create this: int* a = new int[500000000]; however anything larger then std::array<int, 2000000> a; causes a segemtation fault. Hence I am not sure about the properties of the std::array – Vincent Oct 31 '14 at 17:24
• But of course the both of you are right in addressing the issue. – Vincent Oct 31 '14 at 17:28
• @Vincent Blowing the stack is indeed a problem with stack-allocated memory. That's why linear algebra libraries often provide static and dynaic matrix classes. One other solution would be to stack-allocate small matrices and the heap-allocate the big ones; the size of the matrix would be used to chose the allocation policy :) – Morwenn Nov 1 '14 at 10:28
• @Morwenn: excellent suggestion. Well it seems that the more you think about the details of the implementation the actual class becomes more and more complex. In fact in response to your suggestions I posted this stackoverflow.com/questions/26680786/… question yesterday and it seems that there are actually quite a few options concerning that allocation of the stack/heap size. Personally I have not tried out the performance differences yet, but I certainly will. – Vincent Nov 1 '14 at 10:59