 # Proposal: meta-programmed matrix operations?

Hi,

I recently needed a few matrix operations for a filter where I am using ranges. In order to be able to feed those objects in `vtkMath`, there needs to be templated versions of each operation. It is pretty straightforward to do so with most methods, like `vtkMath::Dot` for instance.

I think that this would be the occasion to do loop-free meta-programmed matrix operations, like `vtkMath::MultiplyMatrix` for instance, templated on the matrices dimensions and on the parameters as well as the matrix / vector parameters.

Is anyone against such a feature?

Could you please provide examples that show how things are done now and how could be done with the new templates?

Right now, you can do matrix multiplication this way:

``````// A = 1 0
//     0 1
const double A = { { 1, 0 }, { 1, 0} };
// B = 1 0 0
//     0 1 0
const double B = { { 1, 0, 0 }, { 0, 1, 0 } };

double C
vtkMath::MultiplyMatrix(A, B, 2 /* rowA*/, 2 /*colA*/, 2 /*rowB*/, 3 /*colB*/, C)
``````

You need to put double** or float** as inputs and output. Now, using ranges, one could do that:

``````vtkNew<vtkDoubleArray> arrayOfMatricesA, arrayOfMatricesB, arrayOfMatricesC;
/**
* Depending on the matrix dimensions, fill the arrays.
* One could fill matrices row by row.
*/
// rangeA is a range of 2x2 matrices
auto rangeA = vtk::DataArrayTupleRange<4>(arrayOfMatricesA)
// rangeB is a range of 2x3 matrices
auto rangeB = vtk::DataArrayTupleRange<6>(arrayOfMatricesB);
// rangeC is a range of 2x3 matrices
auto rangeC = vtk::DataArrayTupleRange<6>(arrayOfMatricesC);
using Scalar = typename decltype<rangeC>::value_type;

for (vtkIdType id = 0; id < arrayOfMatricesA->GetNumberOfTuples(); ++id)
{
// transposeA and transposeB would be false by default, but I put them for completeness
vtkMath::MultiplyMatrix<Scalar, 2 /*rowA*/, 2 /*colA == rowB*/, 3 /*colB*/,
false /*transposedA*/, false /*transposedB*/>
(rangeA[id], rangeB[id], rangeC[id]);
}
``````

You could use the first example as well by feeding `*A`, `*B`, and `*C` in this templated version.

I just realized that `Scalar` can be induced with meta-programming, so I think it would not be necessary.

I implemented matrix multiplications in `vtkMath` in master. Look for `MultiplyMatrix`, `MultiplyMatrixWithVector`, and `Dot`.

Matrices are 1D vector that get rearranged in row-column fashion depending on the template parameters.

You can change the layout of matrices at compile-time using some template flags. 3 layouts are supported:

• vtkMath::MatrixLayout::Identity - matrix is read as is
• vtkMath::MatrixLayout::Transpose - the transposed matrix is read when multiplying
• vtkMath::MatrixLayout::Diag - the matrix is a diagonal matrix, i.e. values of the input array are values of the diagonal.