Implementation of the basic matrix operations for embedded application in C++

Question

I have been developing a control software in C++ and for implementation of the control algorithms I need basic matrix operations like addition, subtraction, multiplication and multiplication by scalar.

Due to the fact that the application is the real time application I need to avoid the dynamic memory allocation so I have decided to exploit the C++ templated class and the nontype parameters

template <class T, uint8_t ROWS, uint8_t COLUMNS>
class Matrix
{
public:
  T array[ROWS][COLUMNS];

  Matrix<T, ROWS, COLUMNS> operator+(const Matrix<T, ROWS, COLUMNS> &m) const
  {
    Matrix<T, ROWS, COLUMNS> result;

    for (uint8_t row = 0; row < ROWS; row++) {
      for (uint8_t column = 0; column < COLUMNS; column++) {
        result.array[row][column] = array[row][column] + m.array[row][column];
      }
    }

    return result;
  }

  Matrix<T, ROWS, COLUMNS> operator-(const Matrix<T, ROWS, COLUMNS> &m) const
  {
    Matrix<T, ROWS, COLUMNS> result;

    for (uint8_t row = 0; row < ROWS; row++) {
      for (uint8_t column = 0; column < COLUMNS; column++) {
        result.array[row][column] = array[row][column] - m.array[row][column];
      }
    }

    return result;
  }

  template <uint8_t N>
  Matrix<T, ROWS, N> operator*(const Matrix<T, COLUMNS, N> &m) const
  {
    Matrix<T, ROWS, N> result;

    for (uint8_t row = 0; row < ROWS; row++) {
      for (uint8_t column = 0; column < N; column++) {
        result.array[row][column] = 0;
        for (uint8_t element = 0; element < COLUMNS; element++) {
          result.array[row][column] +=
              array[row][element] * m.array[element][column];
        }
      }
    }

    return result;
  }

  friend Matrix<T, ROWS, COLUMNS> operator*(double k,
                                            const Matrix<T, ROWS, COLUMNS> &m)
  {
    Matrix<T, ROWS, COLUMNS> result;

    for (uint8_t row = 0; row < ROWS; row++) {
      for (uint8_t column = 0; column < COLUMNS; column++) {
        result.array[row][column] = k * m.array[row][column];
      }
    }

    return result;
  }

  friend Matrix<T, ROWS, COLUMNS> operator*(const Matrix<T, ROWS, COLUMNS> &m,
                                            double k)
  {
    return k*m;
  }
};

I have doubts regarding the decision to have the matrix itself i.e. the array as a public member of the Matrix classes.

Answer 1

I have doubts regarding the decision to have the matrix itself i.e. the array as a public member

I don't like having the array as a public member because it means accessing individual elements looks like a bit strange compared to the other operations.

int main()
{
    Matrix<int, 5, 5>   m1;

    // this looks funny.
    std::cout << m1.array[2][2] << "\n";

    // I would want to access the elements like this:
    std::cout << m1[2][2] << "\n";
}

I would add an operator[] to your class to handle the array access directly.

---

I need to avoid the dynamic memory

Not using dynamic memory allocation can quite easily make the application much slower. You need to make sure you are not copying whole arrays around when a dynamically allocated resource can be moved around with a single pointer.

  // This operation returns by value.
  // Which means that unless the compiler finds an optimization
  // You need to copy the result out of the function.
  Matrix<T, ROWS, COLUMNS> operator+(const Matrix<T, ROWS, COLUMNS> &m) const

  // Since you used a std::array this means copying all the
  // individual elements from source to destination.

If you had used std::vector (ie. some dynamic memory allocation) then you can move the content out of the function. This means internally you only have to move a single pointer out of the function. For large arrays the difference can be significant (copying thousands of 32 bit values or copying one 32 bit value).

Now "modern" common C++ compiler have some pretty good optimizations but if you are on an embedded platform these optimizations can be spotty, and you would need to validate how good your compiler is and if it can do RVO or NRVO (Return Value Optimization or Named Return Value Optimization).

---

Side note.

There are some really good optimized matrix libraries out there. Have you tried using one of those.

I am by no means an expert but I have seen the code (a long time ago) where they don't do any of the operations until they actually need the result.

i.e. Matrix Addition/multiplication etc build an expression tree. The expression tree is only evaluated when you try and access an element of the result. This means you can do optimizations on the full expression (like multiplication by an identity can be removed, multiplication by zero can drop a whole sub-tree of the expression etc). This means you only do the operations that actually get you the values you need and don't brute force do all operations.