C++ Simple matrix for the rest of us

The Introduction

“…There are issues such as dynamic memory allocation in C/C++, which is forbidden, under the DO-178B standard, in safety-critical embedded avionics code…”
The source

Before knowing the above, the other day I have developed a handy little template to hold and provide, a plain 2d array, completely on the stack (aka the “matrix”).

It is also fully “compile-time capable” with the full static internal presentation.

It is very handy to gently enforce company-wide C++ policies and to give some useful functionality in return. I am calling this approach “almost a pod”.

Before we proceed: I do know about Expression Templates and few other optimization techniques. Alas. I do consider them an overkill in this context.

All of the dbj::stack_matrix is in this one header file.

Tests are in this another, single header file.

Probably the main reason for having matrices in some app is matrix multiplication. With this API it is rather silly simple, one could say:

taken from

https://github.com/dbj-systems/dbjplusplus/blob/master/test/dbj_static_matrix_test.h

/* make 3 matrices */

using A = stack_matrix<int, 3, 3, DBJ_UID >;

using B = stack_matrix<int, A::cols(),3,DBJ_UID >;

using R = stack_matrix<int, A::rows(), B::cols(), DBJ_UID >;

/* fill them with a callback */

A::for_each(filler);

B::for_each(filler);

/* multiply them*/

stack_matrix_multiply<A, B, R>();

R contains the result

After multiplication, A,B and R relationship is: A[n][m] x B[m][p] = R[n][p] . Columns of A must be the same as rows of B, and R matrix size must be A rows x B cols.

Nice, simple, functional and easy. Since this is all stack solution, the maximum matrix size matter here a lot and is approx 64KB.

This API is more than enough for the 99% of real-life use cases. No huge and complex matrix libraries required.

In case you need matrix addition, difference and such I am sure you will find it very easy to add them by observing how is multiplication implemented, in there.

This is, after all, an open-source. Just please respect the copyright.

In case you just want to use ‘dbj::static_matrix’ as matrix storage, you are most welcome. You will be also more than capable to do so. By looking into the two headers, whose GitHub links are provided.

Example. I personally use its storage in a C++ standard way

// unique int[5][7] on the stack

using A = stack_matrix<int, 5, 7, DBJ_UID >;

// make the standard C++ reference wrapper

std::reference_wrapper<int[A::rows()][A::cols()]>

ref_a = std::ref(A::data());

// or make a native reference to the matrix inside

// int(&)[5][7]

A::matrix_ref_type mr_a = ref_a;

That way I can reach the storage inside and use it elsewhere, to change it or to read it.

If I am using some external API written in a proper standard C++

// external API

template<

typename T, size_t N, size_t M, size_t P

void multiply(

T(&a)[N][M],

T(&b)[M][P],

T(&c)[N][P]

);

I can use this in a “super silly” way:

1 2	// use the 3rd party API and dbj static_matrix storage multiply(A::data(), B::data(), R::data());

Really, a celebration of modern C++. But.

[nextpage title=”Why yet another Matrix?”]
Here is why I did make it.

In standard C++ community 2d array’s is the source of endless ongoing debate. Especially the large ones created on the heap. While at the same time, approx 90% of use cases, in standard programs, do require, a plain old 2d array

// C/C++ plain old data structure
int matrix[3][3];

But, It is a bit dangerous to let your team use this “naked”. It is not that comfortable if you use them regularly and it is dangerous since the stack space can deplete very fast.

As visible in the code behind this API, standard C++ is giving us quite a few mechanisms and features to encapsulate the rules of using such simple native matrices. Making it

Fast Safe and Comfortable

Actual usage might be a bit unusual because instances are not necessary. Just the types. Example

// require matrix as a type

// the rest as usual

template< typename MX, typename T >

void test_mx ( T new_val)

{

// printarr is a method on

// the stack_matrix type

MX::printarr(your_print_function);

}

To use the test function declared as above, we pass our matrix as a type and the rest of arguments as value:

constexpr size_t R = 3, C = 3;

using mx9a

= stack_matrix<int, R, C, DBJ_UID >;

using mx9b

= stack_matrix<int, R, C, DBJ_UID >;

// make sure they are different

// at compile time

static_assert(

mx9a::uuid() != mx9b::uuid()

);

// us the function above

// pass the matrix

// as template argument

test_mx<mx9a>( 0 );

test_mx<mx9b>( 100 );

This is a very fast solution

Above means, there are no copy or move mechanisms being used whatsoever, No complex, computation-intensive, matrix copying or swapping. The matrix simply stays in one place, where it was first created. Much like std::array conceptually works.

If one needs, to pass instances of this matrix in and out of function, as values, this is perfectly unnecessary but perfectly doable:

// note: we receive argument as value

// we return the matrix as value

// none of it does matter.

// There are no instance

// data inside, compiler has no job

// to try copy/move elision

// there is nothing to copy or move

template< typename MX>

MX test_mx_arg_retval(MX the_mx)

{

// change a cell

the_mx.data(

the_mx.rows() - 1,

the_mx.cols() - 1) = 1234;

// for compiler this is

// "zero effort" operation

return the_mx;

}

And the usage of the above function follows:

// yes static_matrix has a default constructor

// which does nothing by the way

mx9a mxa = test_mx_arg_retval(mx9a());

mx9b mxb = test_mx_arg_retval(mx9b());

mxa.printarr(your_print_function);

mxb.printarr(your_print_function);

We receive argument as value, we return the matrix as value. It totally does not matter. There is no instance bound data inside, Compiler has no job to do. Copy/move elision is not applicable. There is nothing to copy or move. No instances are necessary, get used to it :)

There are more usage samples, all in the testing code provided.

The ‘your_print_function’ from above (same as dbj::console::print) is this

inline auto your_print_function

= []( auto const & first_param, auto const & ... params)

{

std::cout << first_param ;

// if there are more params

if constexpr (sizeof...(params) > 0) {

// recurse

print(params...);

}

return print;

};

[nextpage title=”Are we done here, then?”]

No. We are not done here yet. I would like to turn your attention to the problem hidden inside and its solution in dbj::static_matrix .

A general problem description.

If everything is static inside the type it is all shared, between all instances of the same type. Consider this:

// this is not a type

// this is a template

template<typename C> struct buffer {

// every instance of template

// will share this value

static inline mutable size_t count = 1024U ;

};

If we create a few template definitions we get the same number of unique types, like so:

// 3 concrete unique types

using c_buffer = buffer<char>;

using w_buffer = buffer<wchar>;

using n_buffer = buffer<int> ;

So far so good. But. count is shared for all the instances of the same type.

// two instance of the same type

c_buffer instance_1 ;

c_buffer instance_2 ;

instance_1.count = 42 ;

// true

bool the_same_value = instance_2.count == instance_1.count ;

Since the count is shared on the type level (it is static type member) one change applies to every other instance of that type. That is not so good unless it is a design requirement.

Same problem in the context of dbj::static_matrix

Let us assume dbj::static_matrix has this problem. Let us assume that is the issue:

// the unique type

// but made from the wrong stack matrix

using mx9 = stack_matrix<int, R, C>;

// two instance of the same type

// made on the heap

mx9 * mt1 = new mx9;

mx9 * mt2 = new mx9;

// get the reference to cell [1][1] from the type

// because we can

mx9::matrix_ref_type mx = mx9::data();

// give it value 1

mx[1][1] = 1;

// print all three of them

mx9::printarr( dbj::console::print );

mt1->printarr( dbj::console::print );

mt2->printarr( dbj::console::print );

// they are all the same ...

assert(mx[1][1] == mt1->data()[1][1] == mt2->()[1][1])

We are getting 3 identical matrix outputs. This is because there is no instance bound data inside, it is all shared between all instance of the unique type mx9 There is a single data member inside and it is static, so-called “class-wide”

// just a pod 2D array on the stack
inline static value_type data_[R][C]{};

Thus each instance of the template definition (aka the type) for the same template arguments T,R,Cthe combination will share the same 2d array hidden inside.

What matters here are no objects (aka instances) but types aka classes. Which are in this case template instantiations.

Differentiate between same types

What are we going to do? We need, for example, two 3×3 matrices of int’s but we could not rely on instances to keep them separate in the memory.

We could redesign and re-do this template so that nothing is static, true. But this will lead to completely different run-time semantics. And we want, compile-time. We want our matrices as fast as possible, as simple as possible, and not made on the (slower) heap. Here is the solution.

Make different template definitions

Add type Id as the template argument, and use it to make them definitions different.

Thus, template arguments, on our matrix, will be these:

template <

typename T,

size_t R,

size_t C,

unsigned long UID_

class stack_matrix final ;

There is a new template argument:UID_ It is quite enough for the UID_ to be just an ever increasing number. When we create template definitions (aka types), simplified view on what we have now is this:

// two differnt types, thanks

// to the UID template argument

using mx9a =

stack_matrix<int, 3, 3, 1 >;

using mx9b =

stack_matrix<int, 3, 3, 2 >;

The last template argument effectively creates two distinct types above, which are actually both the same ‘thing’: A static matrix, of size 3×3 with integers in the cells. On the stack, all static. Simple, fast and nothing on the heap. But two different types.

The implementation

The UID_ argument. This is what makes it possible. And, UID_ has to be a compile-time constant. My quick (but not dirty) implementation:

#define DBJ_UID __COUNTER__ + 100

Where __COUNTER__ is MSVC/GCC/CLANG predefined macro value giving 101,102,103 … on the level of the compilation unit.

// dbj stack static matrix solution

// two idfferent types

// separated by the UID_ value argument

using mx9a =

stack_matrix<int, R, C, DBJ_UID >;

using mx9b =

stack_matrix<int, R, C, DBJ_UID >;

// different

static_assert(mx9a::uuid() != mx9b::uuid());

Application-wide solution

I assume we all know what the “compilation unit” is in C/C++. Sometimes called “translation unit”. I also assume you do the following every time.

Use the “all is in the headers” approach.In essence modern variant of “Single Compilation Unit” application composition method. In this scenario, the__COUNTER__ macro will produce ID’s unique on the level of application. Thus all the types of the dbj::static_matrix<T,R,C, DBJ_UID> will contain the unique ID’s.
Use namespaces.In this scenario even if you do not use technique 1 you will have separate types because they will be adorned with namespaces.static_assert( ! std::is_same_v <space_a::mx9 , space_b::mx9 >);

The two above or just the second technique will solve all the problems on the application level. In the very unlikely event, you have them.

Conclusion

As ever if anybody needs advice on how to use dbj matrix, in her project please do let us know, in the comments below, or through any other modern communication channel available.

Enjoy the standard C++