Comment by chancho | Hacker News Reader

Comment by chancho | Hacker News Reader

Parent Full thread

Fortran's (alleged) dominance in scientific is probably attributable to tradition (they still teach it to undergrads in non-CS departments) but the native multidimensional arrays have a much bigger impact than aliasing. In Fortran you just index your array like A(i,j,k) and the compiler will compute (and optimize) the addressing for you. In C, a typical (non-computer) scientist who doesn't really focus on mundane shit like this will end up writing something like

    for ( int i=0; i<ni; ++i )
    for ( int j=0; j<nj; ++j )
    for ( int k=0; k<nk; ++k ) {
       a[ k*ni*nj + j*ni + i ] = ...
    }

which sucks. Optimizing multidimensional array access (which characterizes most of scientific computing) is much easier for a Fortran compiler.

This looks like a lot of calcs for the inner loop but is quite easily optimized. The compiler knows that (j * ni + i) is constant in the k loop and that ni * nj is constant. Check the assembly output. But it is probably better to reorder the loops and traverse linearly through memory so that each cache line brought down is fully consumed in order.

The point of that code example is not the calculations, but accessing memory in a cache-friendly way.

Well, yeah, but since Fortran stores things in column-major that's not really an issue. My broader point was that all of these details are hidden from programmers, so there's less for the programmer to screw up and more room for the compiler to work within.

AFAIK, if you want to work with BLAS, LAPACK and friends you have to structure your matrices like this (as a big block of data) rather than pointers to pointers to pointers (ie a[i][j][k]). So even if there is some inefficiency in filling the data structure, presumably the actual matrix multiply or SVD or what-have-you, is actually quite fast.

This is what people mean when they say pointers may be able to look like arrays, but arrays are not necessarily pointers.

Consider, "float a[X][Y][Z]" and "a[i][j][k]".

I don't know how widespread this is, but my uni teaches Python to the physics (and AFAIK most science) undergrads. Some also do C, but no Fortran.

Reply on news.ycombinator.com