Module 10: Matrices

• how do computers store 2D matrix data?
• what is the complexity of multiplying matrices?

Much of the data we take from the world around us stores quite naturally into one-dimensional arrays. For example, computer sound files really just contain arrays of numbers representing the strength of the sound signal at various times. A waveform file contains about 44,000 data points for each second of sound (twice that much if it’s recorded in stereo). MP3 files are similar except that compression techniques are used to make the amount of data smaller.

Other kinds of data, however, require more dimensions. Picture files (gifs, bitmaps or jpegs) all contain data arranged into two-dimensional arrays of pixels. Weather models require three dimensional arrays of cells (2D location altitude), as do many engineering models (stress in a solid object such as an automotive casting).

So, can we do multi-dimensional arrays in C?

The answer is yes and no. C has a way to declare multi-dimensional arrays, but, except in very limited circumstances, there is no useful way to pass multi-dimensional arrays to functions. So, we’ll show you how to do it using single dimensional arrays.

Mapping 2D into 1D

Matrices are two dimensional arrays of numbers. Here are two views of the same 2D matrix, one with colouring to add clarity to the rows:

The matrix has three rows of four columns each. We can number its rows 0 through 2 (counting from zero, as we like to do when manipulating arrays) and its columns 0 through 3. We can refer to any location in the matrix using both the row and the column: in mathematics, we would say that $A_{i,j}$ refers to the value in row $i$ and column $j$ of this matrix. We might like to do similar things in programming, but the computer’s memory doesn’t really work that way: memory is one-dimensional. That is, memory has addresses that range from zero to some large number (like $16 \times 1024^3 - 1$ when a computer has 16 GiB of memory). When we want to represent matrices in memory, we need to find a way to represent a 2D matrix with a 1D array.

You might be thinking, “what about programming language X (e.g., Matlab), which gives me two-dimensional matrices?” The answer is that, while there might be an appearance of 2D matrices, the real computation “under the hood” by any efficient numerical computation software is being done using the kinds of techniques we describe below.

We can build a one dimensional array out of a two dimensional one by simply storing one row after another within an 1D array of length rows $\times$ cols. When rows are stored together in this row, moving to the next row requires jumping ahead by cols spaces: that’s how wide a single row is. When moving to the next column, we only need to move by one space within a row, as shown below:

To find the 1D position of any element $(i,j)$ from the 2D matrix, we can calculate an index as i * cols + j. For example, the value stored at $(0,3)$ in the matrix will be stored at index $0 \times 4 + 3 = 3$ in the 1D array; $(2,2)$ can be found at $2 \times 4 + 2 = 10$. We can also go the other way. Given an index for the 1D array, we can compute (i,i) using:

i = index / cols
j = index % cols


Matrix multiplication

Now let us consider a more complex example—matrix multiplication. If you recall your linear algebra, you remember that matrix-matrix multiplication requires that the number of columns of the first matrix (call it A) must be equal to the number of rows of the second (B). The product matrix P will then have dimensions equal to the number of rows of A by the number of columns of B. The value of each element in the product matrix P is computed by multiplying the rows of A by the columns of B and summing the elements.

In the lecture, we used matrix multiplication as an example of a design problem that required us to follow three steps:

problem analysis
figuring out how the problem works; how we might solve it
algorithm design
figuring out how to break our solution down into simple, repeatable steps; how a computer might solve it
implementation
converting our algorithm (expressed as pseudocode) into a language like C++; how a compiler can read our solution

If you missed this lecture (or you need a refresher), I highly recommend that you watch the lecture capture and then reproduce that process. At the end of the lecture, we had the following results:

Analysis

When adding matrices $A + B = C$, we can calculate the value of any element in the result $C_{i,j}$ as:

$$C_{i,j} = \sum_{k=0}^{n} A_{i,j} \times B_{k,j}$$

Algorithm

for i in [0, rows):
for j in [0, cols):
sum = 0
for k in [0, n):
sum += a[i,k] * b[k,j]
c[i,j] = sum


Implementation

void matrixMultiply(double a[], double b[], double c[], int rows, int cols, int n)
{
int rowlenA = n;
int rowlenB = cols;
int rowlenC = cols;

for (int i = 0; i < rows; i++)
{
for (int j = 0; j < cols; j++)
{
double sum = 0;
for (int k = 0; k < n; k++)
{
sum += a[i*rowlenA + k] * b[k*rowlenB + j];
}
c[i*rowlenC + j] = sum;
}
}
}


(c) 2009â€“2016 Michael Bruce-Lockhart, Theo Norvell, Dennis Peters and Jonathan Anderson. Licensed under a Creative Commons Attributionâ€“Noncommercialâ€“Share-Alike 2.5 Canada License. Permissions beyond the scope of this license may be available at theteachingmachine.org.