Here's my code for matrix multiplication, notice that the size of the matrices is checked compile-time.

This is actually the beginning of a templated matrix class.

Code:

template<int rows, int cols, typename T=double>
class matrix
{
public:
T& operator()(int r, int c)
{
return data[r-1][c-1];
}
private:
T data[rows][cols];
};
template<int rowsA, int colsA, int colsB, typename typeA, typename typeB>
matrix<rowsA,colsB,typeA> operator*( matrix<rowsA,colsA,typeA> A, matrix<colsA,colsB,typeB> B )
{
matrix<rowsA,colsB,typeA> C;
//For each element in C
for (int currentCRow = 1; currentCRow <= rowsA; ++currentCRow)
for (int currentCCol = 1; currentCCol <= colsB; ++currentCCol)
{
typeA sum = 0;
for (int n=1; n <= colsA; ++n)
sum += A(currentCRow,n) * B(n,currentCCol);
C(currentCRow,currentCCol) = sum;
}
return C;
}
template<int rows, int cols, typename T>
std::ostream& operator<<(std::ostream& out, matrix<rows,cols,T> M)
{
for (int r=1;r<=rows;++r)
{
out << "(";
for (int c=1;c<cols;++c)
out << M(r,c) << ",";
out << M(r,cols) << ")\n";
}
return out;
}

Some example code

Code:

int main()
{
matrix<2,3> A;
matrix<3,2> B;
A(1,1)=1;A(1,2)=2;A(1,3)=3;
A(2,1)=4;A(2,2)=5;A(2,3)=6;
B(1,1)=7;B(1,2)=8;
B(2,1)=9;B(2,2)=0;
B(3,1)=1;B(3,2)=2;
// (1 2 3)
// A = (4 5 6)
// (7 8)
// B = (9 0)
// (1 2)
cout << A*B;
}