Row-major or column major?

[Chuck_Hans]

It is common sense that FORTRAN store array in column-major order. But many functions, libraries treat array in row-major order. For example, in ifort document for matmul:

A is matrix

  [ 2  3  4 ]
  [ 3  4  5 ],
B is matrix

  [ 2  3 ]
  [ 3  4 ]
  [ 4  5 ],

The result of MATMUL (A, B) is the matrix-matrix product AB with the value

  [ 29  38 ]
  [ 38  50 ].

Apparently, the matrix multiplication is done in row-major order.

So, if I have to store an array in column-major order, but I do the array operations in row-major order, isn't that confusing? Or perhaps I just store array in row-major order?

[Anonymous66]

Wouldn't it be more confusing if Fortran broke with the standard conventions for matmul?

Row-majororder vs. column-major order is simplya question of how elements of an array are stored. It has no effect how an array index's meaning.

[Chuck_Hans]

Yes, it would be more confusing if it is against standard convention.

FORTRAN treat arry in column-major order, so when handlling an array, I would do it like this:
[fortran]do i = 1, m do j = 1, n arr(j, i) = ... end do end do[/fortran] in other words, I reverse the row, column index, not just for storage, but also in subsequent calls. Part of the reason for this is for efficiency(although maybe not too much).

As of your reply, do you mean that column-major order is just for storage? When treating arrays, do it like they are in row-major order.

[mecej4]

> Apparently, the matrix multiplication is done in row-major order

The notion of a matrix, and the mathematical properties of matrices, were well defined and widely used long before the transistor was invented. These notions have not changed after the advent of the computer. Therefore, the computer language designer has no say in defining matrix multiplication. Pictorial representations such as rows, columns and tables, while popular and helpful, are not essential in defining abstract mathematical concepts.

Row-major and Column-major are conventions for mapping a two-dimensional array to a one-dimensional array; that is, translating a matrix name and two indices to a memory address.

[Anonymous66]

Row-major order vs. column-major makes a difference for traversing a matrix because for largearrays it is more efficientto access elements ofthearray in the order they arestored.It is only a question of what order the elements of an array are accessed, not of the array's shape or the index's meaning.

If arr is ann by m matrix, to store values to it element wise, write:

[fortran]do j = 1, n do i = 1, m arr(i, j) = ... end do end do[/fortran]
The order of the nested do loops is the only thing that changes when dealing with arrays stored in row-major order vs. column major order.

[Chuck_Hans]

In fact, you are making me more confusing.

The example you give is an m by n matrix, not an n by m matrix.

Like you said, "for largearrays it is more efficientto access elements ofthearray in the order they arestored". If arrays are stored in column-major order, for function matmul, despite convension, if it is in column-major order, would it be more efficient?

For a simple case, there are two m by n by n arrays, mat_a(m, n, n), mat_b(m, n, n). The operation to do is an matrix multiplication on mat_a(2,:,:) and mat_b(3,:,:). In row-major scheme, it would be:
[fortran]Integer :: mat_a(m, n, n), mat_b(m, n, n), mat_c(n, n) Integer :: i, j, k Do i=1,n Do j=1,n mat_c(i,j) = 0.d0 Do k=1,n mat_c(i,j)=mat_c(i,j)+mat_a(2,i,k)*mat_b(3,k,j) End Do End Do End Do [/fortran] In columb-major order, it would be
[fortran]Integer :: mat_a(n, n, m), mat_b(n, n, m), mat_c(n, n) Integer :: i, j, k Do i=1,n Do j=1,n mat_c(i,j) = 0.d0 Do k=1,n mat_c(i,j)=mat_c(i,j)+mat_a(k,i,2)*mat_b(j,k,3) End Do End Do End Do [/fortran]Or if we store dimension m in column-major order, n by n in row-major order:
[fortran]Integer :: mat_a(n, n, m), mat_b(n, n, m), mat_c(n, n) Integer :: i, j, k c = matmul(mat_a(1,1,2), mat_mul(1,1,3))[/fortran]Which scheme should I take? As a matter of fact, this situation is what I encounter now. This is why I am confused.

[mecej4]

One reason for the confusion may be your choice of matrix multiplication as the basis for the discussion of column-major versus row-major storage. Matrix multiplication is rather peculiar and poses the following obstacle to achieving efficienty. Each element of the product matrix C = A x B is formed as the inner product of a row of A with a column of B. With multitiered memory (i.e., with cache(s)), this causes problems because one of the vectors appearing in the inner product is going to be accessed with non-unit-stride, no matter which of the two storage schemes is used.

If column-major storage is used, during the computation of A^TB both vectors are accessed with unit-stride. Likewise, if row-major storage is used, AB^T involves unit-stride access. These combinations, are, however not as frequent as the combination AB.

Decades ago, memory was very slow and was sometimes arranged in banks, so the effect of stride was somewhat different.

[styc]

What causes the awkwardness is the interpretation of matrix multiplication as dot products. If it is interpreted as weight sums of the columns of A with weights taken from columns of B, then no strided accesses will show up.

[mecej4]

That is a valid criticism of what I wrote, having adopted the viewpoint of the originator of this thread, that is, "compute one element of the product at a time".

Different views of matrix multiplication are covered in the book Matrix Computations by Golub and van Loan.