I implemented a fixed-size matrix class which supports basic matrix operations (cannot use C++11). What do you think, and how could it be improved?#ifndef FIXED_SIZE_MATRIX_H__#define FIXED_SIZE_MATRIX_H__#include <iostream>#include <cstdlib>template <typename Ty, int N, int M = N>struct FixedSizeMatrix {  typedef Ty value_type;  union {    struct {      Ty element[N][M];    };    struct {      Ty flatten[N * M];    };  };  // Access with bounds checking  Ty& operator()(int r, int c) {    assert(r >= 0 && r < N);    assert(c >= 0 && c < M);    return element[r][c];  }  // Return matrix transpose  FixedSizeMatrix<Ty, M, N> T() const {    FixedSizeMatrix<Ty, M, N> result;    for (int r = 0; r < N; ++r) {      for (int c = 0; c < M; ++c) {        result.element[c][r] = element[r][c];      }    }    return result;  }  // Return matrix inverse  FixedSizeMatrix<Ty, N, M> I();};///// Matrix Inverse///// A design choice was made to keep FixedSizeMatrix an aggregate class to enable// FixedSizeMatrix<float, 2, 2> = {1.0f, 1.0f, 0.5f, 0.2f} initialization.// With c++11 it would be possible to create a FixedSizeMatrix base class// and derive all variations from it, while retaining the brace// initialization.//// Matrix inverse helpersnamespace detail {// General version not implementedtemplate <typename Ty, int N, int M>struct inverse;// Matrix inversion for 2x2 matrixtemplate <typename Ty>struct inverse<Ty, 2, 2> {  FixedSizeMatrix<Ty, 2, 2> operator()(FixedSizeMatrix<Ty, 2, 2> a) {    FixedSizeMatrix<Ty, 2, 2> result;    Ty det =        a.element[0][0] * a.element[1][1] - a.element[0][1] * a.element[1][0];    assert(det != 0);    result.element[0][0] = a.element[1][1] / det;    result.element[1][1] = a.element[0][0] / det;    result.element[0][1] = -a.element[0][1] / det;    result.element[1][0] = -a.element[1][0] / det;    return result;  }};}  // detail// Define matrix inversetemplate <typename Ty, int N, int M>FixedSizeMatrix<Ty, N, M> FixedSizeMatrix<Ty, N, M>::I() {  return detail::inverse<Ty, N, M>()(*this);}/////  Matrix operations///// Matrix producttemplate <typename Ty, int N, int M, int P>FixedSizeMatrix<Ty, N, P> operator*(FixedSizeMatrix<Ty, N, M> a,                                    FixedSizeMatrix<Ty, M, P> b) {  FixedSizeMatrix<Ty, N, P> result;  for (int r = 0; r < N; ++r) {    for (int c = 0; c < P; ++c) {      Ty accum = Ty(0);      for (int i = 0; i < M; ++i) {        accum += a.element[r][i] * b.element[i][c];      }      result.element[r][c] = accum;    }  }  return result;}// Unary negationtemplate <typename Ty, int N, int M>FixedSizeMatrix<Ty, N, M> operator-(FixedSizeMatrix<Ty, N, M> a) {  FixedSizeMatrix<Ty, N, M> result;  for (int e = 0; e < N * M; ++e) result.flatten[e] = -a.flatten[e];  return result;}#define MATRIX_WITH_MATRIX_OPERATOR(op_symbol, op)                            \  template <typename Ty, int N, int M>                                        \  FixedSizeMatrix<Ty, N, M> operator op_symbol(FixedSizeMatrix<Ty, N, M> a,   \                                               FixedSizeMatrix<Ty, N, M> b) { \    FixedSizeMatrix<Ty, N, M> result;                                         \    for (int e = 0; e < N * M; ++e)                                           \      result.flatten[e] = a.flatten[e] op b.flatten[e];                       \    return result;                                                            \  }MATRIX_WITH_MATRIX_OPERATOR(+, +);  MATRIX_WITH_MATRIX_OPERATOR(-, -); #undef MATRIX_WITH_MATRIX_OPERATOR#define MATRIX_WITH_SCALAR_OPERATOR(op_symbol, op)                          \  template <typename Ty, int N, int M>                                      \  FixedSizeMatrix<Ty, N, M> operator op_symbol(FixedSizeMatrix<Ty, N, M> a, \                                               Ty scalar) {                 \    FixedSizeMatrix<Ty, N, M> result;                                       \    for (int e = 0; e < N * M; ++e)                                         \        result.flatten[e] = a.flatten[e] op scalar;                         \    return result;                                                          \  }MATRIX_WITH_SCALAR_OPERATOR(+, +);MATRIX_WITH_SCALAR_OPERATOR(-, -);MATRIX_WITH_SCALAR_OPERATOR(*, *);MATRIX_WITH_SCALAR_OPERATOR(/, / );#undef MATRIX_WITH_SCALAR_OPERATORtemplate <typename Ty, int N, int M>FixedSizeMatrix<Ty, N, M> operator+(Ty scalar, FixedSizeMatrix<Ty, N, M> a) {  return a + scalar;}template <typename Ty, int N, int M>FixedSizeMatrix<Ty, N, M> operator*(Ty scalar, FixedSizeMatrix<Ty, N, M> a) {  return a * scalar;}template <typename Ty, int N, int M>FixedSizeMatrix<Ty, N, M> operator-(Ty scalar, FixedSizeMatrix<Ty, N, M> a) {  return -a + scalar;}template <typename Ty, int N>FixedSizeMatrix<Ty, N, N> identity_matrix() {  FixedSizeMatrix<Ty, N, N> result = FixedSizeMatrix<Ty, N, N>();    for (int i = 0; i < N; ++i) result.element[i][i] = Ty(1);  return result;}template <typename Ty, int N, int M>std::ostream& operator<<(std::ostream& out, FixedSizeMatrix<Ty, N, M> a) {  for (int r = 0; r < N; ++r) {    for (int c = 0; c < M; ++c) {      out << a.element[r][c] <<  ;    }    out << std::endl;  }  return out;}#endif
Fixed size matrix implementation
c++;matrix;c++03
Looks good over all, but I have a few thoughts. A few are stylistic, some are technical, and some are just trade offs I want to make sure you're aware of.I would use the value_type typedef for method return types, temporary types, etc. It should make for a bit prettier signature in IDEs and any kind of auto generated documentation.I would try to be consistent with variable names with regards to N/M/r/c. I either use N/M/n/m or R/C/r/c. Some people like to avoid using different capitalizations though, so that might be a valid concern.Template type parameters default to the current ones inside of a templated class, so all of the FixedSizeMatrix<Ty, M, N> types can be simplified to FixedSizeMatrix.It's good that you're consistent, but 2 spaces is rather unusual for C++. I would use either 4 spaces or tabs.You should also have a const version of operator().If you ever decide to make the class more sophisticated, you might regret exposing element. Perhaps hide the struct as a private member, have operator() be unchecked, and have a checked version called at. You could even have a 1 parameter version of each one that would provide the flattened functionality.This unfortunately would break aggregate initialization though, which I suspect is why you've made the decision to have it public.Still though, even if you make keep it public, you can get in a habit of using operator() and at() and then if you ever use C++11, you can properly encapsulate it with minimal code breakage.I would throw a std::out_of_range instead of using asserts in operator(). Note that this and having the checked version named at while using the operator as a pass through is consistent with the behavior of std::vector.Same thing for your other asserts. Exceptions can be caught and handled if desired. Asserts not so much. if you're using asserts for performance reasons, please consider if the tiny difference really matters. Modern C++ compilers have a zero cost exception model where exceptions are only costly if they're actually thrown.Don't use double underscores. You technically can in certain situations, but the rules are very difficult to remember. Your FIXED_SIZE_MATRIX_H__ is technically invalid.Is I common notation for a matrix inverse (makes me think identity matrices)? Even if it is, I would name the method invert. Same thing with the transpose.Any time you're passing the matrix as an argument to a function and you're not modifying it, pass it by constant reference. There's no need to make a potentially very expensive copy just to read from it. There's also the option of having the argument copy but then using the copy to store the result (though that requires that you can operate on the matrix without a dependency cycle between elements). This would allow move semantics in C++11 on temporaries passed as the argument.Everywhere you do Ty(0), you force an interface constraint on Ty in that it must have a constructor that is compatible with a single integer parameter. If you default construct the Ty, you avoid this interface constraint, though you of course add a new one. Consider the two cases this can affect: is a type used more likely to have a default constructor that represents 0 or a single integer constructor? I would lean towards the first one, oddly enough. For example, a Complex might not have an integer constructor.It might be nice to put this in a namespace if you plan on the usage of this to grow to more than just 1 or 2 simple programs. That helps mentally organize your library as one coherent part, and it provides a few technical advantages. Namely, your detail namespace could conflict with someone else's ::detail namespace, and it would be nice to give your macros, though short lived, a prefix of some kind. It's relatively unlikely that anything in this would actually conflict, but if someone were using this with another matrix library, you never know. Namespaces don't cost anything to you as the developer (or to consumers of the class), but they offer a few advantages.