don’t be afraid of variably modified types

This post is about a subject that easily triggers a lot of polemic and diverging opinions, so I will try to be careful to get things right. In any case the message of this post may be summarized as:

• Plain variable length arrays (VLA) might be dangerous for your health, since they may cause stack overflows, and generally are not so easy to handle
• In constrast to that pointers to VLA are a quite different story and may ease the programming with multi-dimensional arrays substantially.

Ancient C only had arrays with a length that was known at compile time. For such arrays this knowledge is deduced by either providing an explicit integer constant expression or by giving an initializer:

double A[3] = { 0 };
double B[] = { 0.0, 0.0, 0.0 };

Here A and B have the same type and contents.

C99 changed this, in that in function scope it now allows to have an integer expression that specifies the length of the array.

enum { L = 3 };
double A[L] = { 0 };              // valid, L is a constant integer expression
double B[] = { 0.0, 0.0, 0.0 };
size_t const K = argc;
double C[K];                       // valid, K is an integer expression, C is VLA
double D[K] = { 0.0, 0.0, 0.0 };  // invalid, D is VLA and can't be initialized 

An array such as C in this example is called a variable length array, VLA. The type of such an array is a variably modified type, VM. We may rewrite part of this example as:

enum { L = 3 };
typedef double fixedArray[L];
fixedArray A = { 0 };
size_t const K = argc;
typedef double variableArray[K];
variableArray C;
variableArray D = { 0.0, 0.0, 0.0 };  // still invalid

Here now variableArray is a VM, and expressions such as sizeof(variableArray) or sizeof C can not be evaluated at compile time but may give different results at each call.

Digression: VLA are considered evil by many programmers and programming language designers. The main argument I see repeatedly is that they may cause stack overflow, SO, in case of errors or used with a large size. A SO is a run time error that is very difficult to diagnose and to debug:

• A function when it is called usually receives its stack location through a special hardware register, the stack pointer, SP. At the beginning of its operation it usually then reserves the space that it needs for its local variable and sets SP to a new value, taking this situation into account. At the end, before returning, the old value of SP is restored.
• If we have VLA, things become complicated. Since the size of a VLA is not known in advance, SP must be adjusted with a variable difference. Since this difference is variable, restoring the old value at the end is not as simple as before: we have to keep track of the amount by which we modified SP. To keep track of that we’d have to place that value on the stack… As you can see not an easy problem to handle, and in particular many things can easily go wrong at run time. A program that has an erroneous SP is one of the worst things that can happen to it. It will write data in places where it shouldn’t it will loose the return address of a function, and so on.

An important feature of VLA is that the variability can be used in all dimensions of a multidimensional array:

double A[n1][n2];
double B[n2][n3];
double C[n1][n3];
.
/* check that the sizes of the matrices correspond */
assert(sizeof(A)/sizeof(A[0]) == sizeof(sizeof(C)/sizeof(C[0])));
assert(sizeof(A[0])/sizeof(A[0][0]) == sizeof(sizeof(B)/sizeof(B[0])));
assert(sizeof(B[0])/sizeof(B[0][0]) == sizeof(sizeof(C[0])/sizeof(C[0][0])));
/* do the multiplication */
for (size_t i1 = 0; i1 < n1) {
  for (size_t i3 = 0; i3 < n3) {
    C[i1][i3] = 0.0;
    for (size_t i2 = 0; i2 < n2)  C[i1][i3] += A[i1][i2] * B[i2][i3];
  }
}

This is a valid code for matrix multiplication, even if none of the bounds n1, n2 or n3 are known at compile time. (The sizeof expressions are ugly, but see below.)

With historical C, such a matrix multiplication would have been more difficult to implement:

• We’d have to have to define the matrices as pointers to double where semantically they are something different for us.
• We’d have to malloc space for them and free this space at the end.
• Without any type safety as guaranteed by the assertions in the above code.
• We’d have to implement the 2-dimensional access operation by ourselves with something like i1 * n1 + n3.

If we don’t want to allocate arbitrary large VLA on the stack, will not be able to avoid the second point, the allocation on the heap. But all other advantages can be kept when using pointers to VLA instead of VLA themselves.

double (*A)[n1][n2] = malloc(sizeof(double[n1][n2]));
double (*B)[n2][n3] = malloc(sizeof(double[n2][n3]));
double (*C)[n1][n3] = malloc(sizeof(double[n1][n3]));
.
for (size_t i1 = 0; i1 < n1) {
  for (size_t i3 = 0; i3 < n3) {
    (*C)[i1][i3] = 0.0;
    for (size_t i2 = 0; i2 < n2)  (*C)[i1][i3] += (*A)[i1][i2] * (*B)[i2][i3];
  }
}
.
free(A); free(B); free(C);

We might perhaps get yourself used to expressions such as (*A)[i1][i2], but if we keep clean spaces around arithmetic operators and assignments, this will not be too difficult I think.

In the first version I insisted to do the computation of the array dimensions with clumsy sizeof expressions. The reason is that the dimensions of or matrices are fixed at the point of their definition. The variables n1 ... may have been changed by some code between the two points, so the only possibility to be sure that we compare the right values is to use sizeof. BTW the same holds for the loop bounds, so ideally we would use sizeof expressions there also.

This is not very beautiful, lengthy and error prone. P99 has the macro P99_ALEN to help us on that. Called as P99_ALEN(A) this does nothing spectacular, it evaluates to sizeof(A)/sizeof(A[0]). But used with an integer N as its second argument it will give you the Nth dimension of the underlying array.

/* check that the sizes of the matrices correspond */
assert(P99_ALEN(A, 1) == P99_ALEN(C, 1));
assert(P99_ALEN(A, 2) == P99_ALEN(B, 1));
assert(P99_ALEN(C, 2) == P99_ALEN(B, 2));
/* do the multiplication */
for (size_t i1 = 0; i1 < P99_ALEN(C, 1)) {
  for (size_t i3 = 0; i3 < P99_ALEN(C, 2)) {
    (*C)[i1][i3] = 0.0;
    for (size_t i2 = 0; i2 < P99_ALEN(A, 2))  (*C)[i1][i3] += (*A)[i1][i2] * (*B)[i2][i3];
  }
}

In a later post we will see how all that can be put together into functions that use VM types easily.