English is currently available.
Fortran DISCUSSION
I initialised a 2 by 3 matrix with a = reshape([1, 2, 3, 4, 5, 6], [2, 3]), expecting the first row to be 1 2 3. Printing row by row gives 1 3 5 and 2 4 6 instead. The same data passed to a C routine also arrives transposed.
Why does Fortran order the elements this way, how do I write the initial values in the row order I see on paper, and does the order of nested do loops over a large array make a measurable difference?
Fortran stores arrays in column-major order: the first index varies fastest in memory. An array constructor is always one-dimensional, and reshape pours its elements into the target shape in that storage order, so 1 and 2 fill column 1, 3 and 4 fill column 2, and 5 and 6 fill column 3. Hence a(1,:) is 1 3 5 and a(2,:) is 2 4 6.
For the same reason print *, a prints 1 2 3 4 5 6, which is memory order, not rows. To print rows, loop over the first index and print a(i, :).
To type the data as it looks on paper, either reshape to the transposed shape and transpose, a = transpose(reshape([1, 2, 3, 4, 5, 6], [3, 2])), or use the order argument: a = reshape([1, 2, 3, 4, 5, 6], [2, 3], order=[2, 1]). Both give rows 1 2 3 and 4 5 6. With order=[2, 1] the second subscript is the one that varies fastest while filling.
For larger tables, continuation lines with one matrix row per source line keep it readable. The square-bracket constructor is Fortran 2003; older code writes (/ 1, 2, 3 /).
Loop order does matter for large arrays. Memory is contiguous along the first index, so the innermost loop should run over the first index: do j = 1, n outside and do i = 1, m inside, touching a(i, j).
A 4000 by 4000 double precision array is 128 MB, larger than the cache of a typical desktop processor. Walking it along columns reads memory sequentially; walking it along rows jumps 32,000 bytes between consecutive accesses and misses the cache almost every time. The slowdown can be several-fold, although an optimising compiler can sometimes interchange simple loops for you. Whole-array expressions such as a = a * 2.0 leave the traversal order to the compiler.
The transposed data in C is the same fact seen from the other side. C is row-major, the last index varies fastest, so a Fortran a(2,3) has the same memory layout as a C a[3][2]: Fortran element a(i,j) is C element a[j-1][i-1]. Nothing needs to be copied; swap the index order in the C code and remember the 1-based versus 0-based indices. NumPy offers the same choice, where arrays created with order='F' use the Fortran layout.
It also affects slicing: a(:, j) is a contiguous column, while a(i, :) is a strided row, which may be copied to a temporary when passed to a routine that expects a contiguous array.
