FORTRAN has multidimensional arrays, and the compilers know about them and optimize aggressively. Most other languages don't. C people just do not get multidimensional arrays. They think an array of arrays is good enough. I tried to convince the Go people to put in multidimensional arrays, with no success.
Years ago, I had an amusing fight with the C++ committee over this.
In C++, you can overload the "[]" operator:
int& operator[] (int i)
But you can't do that for multiple arguments.
int& operator[] (int i, int j)
Why? An obscure feature of C and C++ is that the precedence of the comma is different for function argument lists and subscript lists. In
foo(a,b)
is a function call with two arguments.
foo[a,b]
is an invocation of the comma operator, which ignores the second operand. The syntax charts reflect this. So I proposed that the syntax for comma usage be the same in both
"()" and "[]" forms. If you see "foo[a,b]" in a C/C++ program, it's almost certainly an error. I searched through large volumes of code and was unable to find a single use of the comma operator inside brackets.
But no, there's a "save the comma operator" faction.
> C people just do not get multidimensional arrays. They think an array of arrays is good enough.
I have to say, i dont get what you are trying to say. Do you have a source that explains it?
What is the difference between a multi dimensional array and arrays of arrays? And what does a FORTRAN compiler optimize here? It is a piece of memory written and organized in a certain way. What is there to organize. You basically have a pointer to a memory region with an offset for the position at row X and column Y.
What am I missing? Are you only talking about the syntax for the final address? What has that to do with the compiler? Do you mean the parser?
I'm not a Fortran or C programmer (so take this with a grain of salt), but as I understand it, Fortran is so optimizable because it's arrays are not fundamentally just a memory layout with pointer math. In Fortran, arrays are abstract entities, which allows the compiler to lay them out in the best possible way to take advantage of vectorization for the current code. Originally, pointers weren't even a feature of the language, though they were bolted on in a later version.
Sure, an expert programmer can replicate this speed in C by being careful about memory layouts, etc. But in Fortran, everyone gets it for free.
I would add that since Fortran 90, you can also get scalar-array and array-array operations for free. The `elemental` keyword that automatically broadcasts pure scalar functions to multidimensional arrays works almost in a way reminiscent of NumPy (I'm a weird one who came to do some serious work in Fortran long after known C and Python).
an array of arrays is not necessarily contiguous in C(++). Indeed, if allocating on the heap, you end up with a bunch of discontiguous memory that's not necessarily correctly aligned.
A good tensor implementation accounts for strides that are SIMD compatible (eg; each dimension is a multiple of the SIMD register width).
An array of arrays is necessarily contiguous in C - this is implied by the type. An array of pointers to arrays will, of course, not be contiguous - and is the only way to get a dynamically sized heap-allocated 2D array in C (VLAs give you stack-allocated arrays, with all the size limits that entails).
In C++, this all is best handled by a library class.
Ah, good point. I always forget that VLA types in C99 are actually types, and so you can use them in these contexts as well.
It's a shame they killed VLAs as a mandatory language feature. They didn't make C into Fortran (which I think was the hope, between them, complex numbers, and "restrict"?), but they did make some things a great deal more pleasant to write.
>and is the only way to get a dynamically sized heap-allocated 2D array in C (VLAs give you stack-allocated arrays, with all the size limits that entails).
Why wouldnt you be able to create a dynamic array of arrays with placement new and a cast.
> A good tensor implementation accounts for strides that are SIMD compatible (eg; each dimension is a multiple of the SIMD register width).
Never implemented any tensors, but in my experience, sometimes you better do what GPUs do with some texture formats: switch from linear layout into dense blocks layout. E.g. for dense float32 matrix and SSE code, a good choice is a 2D array of 4x4 blocks: exactly 1 cache line, and only consumes 4 out of 16 registers while being processed i.e. you can do a lot within a block while not hitting any RAM latency.
There's a better solution, to your problem: pass an "index object" as the argument to operator[].
So `T& operator[](index_t)`, where `index_t` could just be a tuple.
a[{1, 2}] = 12;
This has virtually the same syntax, will be compiled in the same way by any modern compiler, and lets you create special index types that have invariants (suppose you want i0 < i1).
I disagree. The tuple is an index into an array. The multiplicity of the tuple indicates (literally, in fact) the multiple dimensions of the array. Ie to the coder, it is the index that takes on the multi-dimensionality, not the array.
It's one of those many tiny things that add up. There's good reasons why C++ gets so much hate and there's good reasons why it's used in so many applications. I see no good reason to not have multiple array indices.
> I see no good reason to not have multiple array indices.
It would add a little bit complexity to the language with only a very smal benefit: The same can be achieved with just a few characters more or () instead of [].
Using Matrix[i,j] is the most intuitive way for most people. It's what you'd use if you don't know anything and it's what you'd understand immediately. That's the benefit for the user.
And to be honest, in the list of things that the C++ folks seems to think is important, multi-index index operators seems like it would rank substantially higher than a lot of the other things in the list in terms of value, and a lot lower in-terms of additional complexity.
If this was true then the BLAS Fortran library would run faster when calling with Fortran as opposed to C, but it doesn't. As someone else pointed out, multi-dimensional arrays are contiguous in C, by definition. It's not like Fortran has some magical powers to unlock instructions only available to that language.
The argument is not that because Fortran has multidimensional arrays it is faster, it is that because Fortran has multidimensional arrays it can and will often be faster in most cases for free.
Obviously an extremely optimised library such as BLAS will run just about as fast in C and Fortran. For a good majority of code out there however people just don't want to put in the effort to carefully optimise array layouts. The feature moves the responsibility of optimisation to the compiler which can in many cases do a better job than a normal programmer with reasonable assurances that the optimisations will not introduce bugs and / or side effects into the code.
Languages like C(++), Rust, and Julia have vector intrinsics that are unfortunately lacking in Fortran. Is there a way to wrap [SLEEF](https://github.com/shibatch/sleef) in Fortran?
While writing fast code in Fortran is easy, I think it is unfortunately harder to write "extremely optimised" libraries like BLAS in Fortran than in these other languages for that reason. Without the low level control, you can't be as explicit about vectorization patterns.
Eg, for BLAS, you normally write a lot of kernels, which you then apply to blocks of the matrices. The blocks will not be contiugous; each column will be separated by the matrix stride (or rows separated by a stride, if row major).
gfortran 8.2 will not succesfully vectorize a matmul kernel (I did file an issue on bugzilla).
C (with vector intrinsics) does great. Julia does well too, but not perfect: there are a couple redundant move instructions per for loop iteration.
When it comes to applying the kernel to blocks, Fortran will also try and make array temporaries of each block, which would cripple performance. In languages with pointers, all you'd have to do is pass a pointer and the stride between columns.
I haven't tried playing with pointers in Fortran, but I guess that would work too?
Julia often refuses to autovectorize anything as soon as you start playing with pointers; does Fortran's stance on aliasing mean it fairs better?
If not, Julia at least has the LLVM vector intrinsics and metaprogramming that let you force it to generate nearly optimal code despite its protests.
But to say something nice about Fortran: dynamically-sized-mutable-stack-allocated-arrays. (-fstack-arrays)
In Julia, stack-allocated objects are immutable, and (partly for that reason) also generally very small. You wouldn't want to create a new 30x30 array every time you want to edit a single value.
You also can't get stack pointers, meaning you can't use masked load/store operations to vectorize code when the array dimensions aren't a multiple of SIMD-vector-width.
This means to write fast code, I normally heap allocate everything. And to avoid triggering the GC, that means keeping the same arrays alive, and avoiding any temporaries.
With Fortran, you don't have to worry about any of that. You can always write convenient and clear code, and it's likely to be very fast. If you hold vectorization in mind while laying out your data and computations, compilers will normally do a great job figuring it out.
Highlights include that that "Stack allocation is a property of the object, mutability is a property of the type", and that mutable objects that do not escape a function are likely to avoid heap allocation (ie, try to avoid calling non-inline functions).
Unfortunately, that can sometimes mean having to roll your own functions like `sum` and `dot_product`. A simple for loop is easy, but that some of these basic functions don't inline by default can make it a little more cumbersome.
Big disclaimer: my experiences are limited, and compiler comparisons are going to be highly variable, depending on the particular code bases.
That said, I have a statistics model that has to be fit thousands (millions?) of times. As a baseline, running it using JAGS (a popular tool for Gibbs sampling) for a given number of iterations and chains in parallel took >160 seconds.
The same model in Julia took 520 ms, g++ & gfortran 480 ms, and ifort took 495 ms.
The code I compiled with g++ used vector intrinsics and SLEEF. Not exactly a revelation that code like that will be fast. But part of my point is that C++ easily lets you take measures like that when you want more performance or control, therefore it is more optimizable. Switching to a commercial compiler wasn't an automatic performance boon.
All else equal, I'd much prefer sticking with an open source compiler suite. I also won't be a student with access to the Intel compilers for much longer.
The actual model is addressing a fun problem (improving the accuracy of satellite-space debris collision probability estimates), and I'm organizing some of it into a "Gibbs Sampling" presentation for tomorrow, so I could put the material online.
Do you mean syntax? Because if layout, then Gorgonia's tensors are actually flat dense layouts controlled by an access pattern. Gorgonia supports fortran layout as well as C layout.
I'm actively working on multiple Fortran projects, including one that's being written from scratch just in the last year in modern Fortran 2018.
Modern Fortran looks completely different from the old school FORTRAN77 most people probably imagine. Gone are the days of fixed-format, where the columns actually mattered, and gone are the days of unreadable ALL CAPITAL LETTERS, and gone are the days of GOTO.
I developed a Fortran implementation of Python's argparse[0] recently to use for another project. The code is nothing like the monstrous spaghetti code of days past---although I do still obsessively line up my code in columns, this isn't important.
In the early 2000's I worked on some Fortran-77 code that we were parallelizing (As a babe, this was my first experience with Fortran-77). We were having strange bugs because it seemed like some MPI calls were just not happening and others were. There were no compiler errors or warnings. Everything was supposedly fine. However, certain properties were propagating across grid elements just fine and others weren't. e.g. A particle might cross a vertical boundary and be just fine in the next CPU's element, but if it crossed a horizontal boundary it would lose it's horizontal momentum but keep everything else.
Eventually we figured out that the lines that weren't working were just over 80 characters long and the ones that were working were under 80 characters long. The compiler we were using was adhering to the Fortran-77 spec of a maximum of 80 characters per line (because that's what would fit on a punch-card), but it was just ripping those lines out and ignoring them without raising any kind of warning. If a line of code was over 80 characters long, it was just silently ignored. All we had to do was split the lines over 80 characters up and that particular bug was squished.
In my humble opinion, if a language can trip you up because of how many characters fit on a punch-card, it qualifies as "old school".
It's basically a switch statement, but without the benefit of descriptive case identifiers. So it's basically,
switch (I-20) {
case 21:
case 22:
...
}
but you have to imagine "21" as the string identifier, not the integer value. `I-20` is indexing the list of labels in the GOTO statement.
The earlier column beginning with
GOTO (1,2,3,4,5,6,7,8,9,10,11,12,13,14,15) I
is easier to read
switch (I) {
case 1:
...
case 2:
...
}
because the integer value of I, the string labels, and their index in the list are all the same.
I've never written Fortran. I did learn to program using TI-85 BASIC, though, and am well versed in C. The GOTO in the Fortran code is possibly the easiest part for me to understand.
If you look at the code, they seem to follow a similar convention. But the computed GOTO described above is using a set of labels employed for what looks like an inner state machine, related to a previous inner state machine (thus the offset by 20). It seems perfectly reasonable (preferable, even) to use an entirely different domain of numeric labels, distinguished by length, for those parts of the code.
What’s the advantage of sticking to the new Fortran? Is it substantially compatible with legacy FORTRAN codebases and systems, despite having significant advances as a language?
Fortran has better aliasing rules than C. Most of the strangeness you have to do with restrict pointers in C are not a concern in high performance code. Because strides and dimensions and types are all laid out explicitly, the optimizer can really go to town. And with co-arrays the same syntax generalizes to multiple compute nodes.
I worked as lead on IBM’s ASTI component (part of xlf) 10 years ago, and kind of miss it sometimes.
Array-based syntax for distributed memory parallel computing (look for coarrays) since Fortran 2008. OOP since Fortran 2003. Whole-array operations since Fortran 90. These are a few of my favorites.
New standards mark obsolescent or deprecate a few ancient features as they come out, however compilers maintain the support for these features. I'd say that the modern Fortran is almost completely backwards compatible with legacy FORTRAN, as far back as FORTRAN IV (predecessor to FORTRAN 77).
> I'd say that the modern Fortran is almost completely backwards compatible with legacy FORTRAN, as far back as FORTRAN IV (predecessor to FORTRAN 77).
That's very impressive! I assumed the compatibility difference was akin to what Perl 6 is to Perl 5- -- that is, not compatible at all.
(1) Compatibility with legacy FORTRAN. Although it's not difficult to work with legacy FORTRAN via C, it's seamless to pass data around between legacy FORTRAN and Modern Fortran.
(2) The mathematical and array syntax. It's tough to beat Fortran's mathematical notation and easy use of multidimensional arrays.
This. As someone who once did a lot of computer graphics programming in C (not C++), I would have killed for proper arrays. Modern C, with variable length arrays, helps, but still...
A long time ago, at a university far, far away, I worked with SGI's gl library, which had FORTRAN bindings. In retrospect, that seems like a very smart idea. At the time, we stuck with the C bindings.
“I don’t know what the programming language of the year 2000 will look like, but I know it will be called FORTRAN.” – Charles Anthony Richard Hoare, circa 1982
It seems every language keeps the same name as it evolves, except Algol/C/C#/Java.
Hoare also said "ALGOL 60 is alive and well and living in FORTRAN 90". I don't have the source of that, but I think it was the answer to someone asking him about it, reported in some Fortran forum.
There's a cute namelist hack that allows parsing arguments as well as serializing them to disk using old, standard language features: https://github.com/RhysU/namelist_cli
I've noticed that Fortran comes up a lot when talking to people who do fluid dynamics simulations.
Curiously, it often does so in context of Python. One specific example I can remember was a library written in Fortran, with unit tests written in Python.
Fortran and Python go together very well (they complement each other in just the right ways and there's all the bindings / data structure compatibility with Numpy that you need). Glue code / UI in Python, Numerics in Numpy + homegrown Fortran, that's how I'd implement a numerical model from scratch today.
Just out of curiosity, what do you use for calling Fortran from Python (e.g. f2py, ctypes)? Do you have any suggestion about how to combine them together (e.g., for parallel calculations)?
I have a friend who writes modern fortran nearly everyday as an actuarial pension analyst in Manhattan. They hired him out of college having never written a "hello world" script. Hope that helps!
Very cool thanks for sharing the project! I used fortran in college for meteorology, had a very oldschool professor. There was another guy who taught Matlab, i liked fortran more honestly. It was the old style fortran , reading thru your repo i think this version does look much more modern.
After writing code in fortran and running on a shared cluster with a buttload of CPU and RAM, working in a software shop doing python science stuff on a mac makes me miss the speed. I didnt give a about how big my source data was , but now with python, occasionally i miss the ol static typed compiled days. I am giddy when i get to do some numpy matrix array type stuff instead of just glueing things together with python magic!
I gather that gfortran can be built on Android, so presumably it's possible. I see conflicting information on whether or not it's actually possible on iOS. I've personally never seen it.
I believe free format came about in Fortran 90, and that FORTRAN 77 was still fixed format--at least that's what a quick Google suggests. But even still, fixed-format F77 has stuck around in lots of legacy code bases for a long time. We only migrated one of our large codebases to free format a few months ago.
I learned Fortran IV and later Fortan 77 in the mid eighties and programmed professionally in Fortran 77 up until 2000. I was taught and programmed with fixed format. I believe that there may have been compiler switches to allow non-fixed format. I'm pretty sure that longer lines were allowed via a compiler switch.
I am still writing new code in fixed format. Because of the wonderful compatibility, I can do it while using a lot of new features from the recent standards which are both supported in GNU Fortran and Intel Fortran. You can call me old school, but I like the constraint of the line length in fixed format, it helps me writing clean code and reduce the "double" loops in my algorithms.
From what I gather, fixed columns was all FORTRAN versions up to F77, since Wikipedia says free format was introduced in Fortran 90 [0]. But it may have been the case perhaps that some compilers supported proprietary directives to support this in earlier versions--I wouldn't know, though.
I understood the Ship of Theseus to have been rebuilt with new materials, not with a new design. I suppose the interpretation is somewhat flexible but I don’t think it quite applies to this case. If a modern program doesn’t compile on an older version, I would interpret that as being a significant change in the language, the referent of “FORTRAN”. I believe GP simply didn’t develop their argument effectively. 1970’s Camaros aren’t interchangeable with 2010’s. To the public they are both Camaros, but to a mechanic they are very different.
>If a modern program doesn’t compile on an older version, I would interpret that as being a significant change in the language, the referent of “FORTRAN”.
A modern C program doesn't compile in old C compilers. That doesn't mean it's not C. You can change tons of ergonomics and keep the language the same. Or you can keep the spirit the same (which in a language are its core semantics and concepts).
>I believe GP simply didn’t develop their argument effectively.
GP just casually said it's a whole new language in the casual, everyday sense, that a e.g. you can say "it's a brand new car now" after you redid the interior and did some repairs on your 10-year old car.
(And I got in with the Ship of Theseus thing to humor the objection -- not that I really consider the GP statement to have to be taken at face value).
Is a 80 year old person the same person and with the same ideas, abilities, look, as when he was a 15 year old person? He has the same name after all.
Same language doesn't mean you can necessarily speak it without learning the new words, forgetting some stuff that went deprecated, learning some syntactic sugar and new features. If you were born in 1400 and lived to today, you'd be going from speaking from Medieval French to modern French with hardly noticing the difference...
Well, I think that's a reasonable analogy, except that with spoken language evolution, the unintelligibility tends to be mutual, whereas with Fortran there's substantial compatibility in one direction (i.e., GFortran will happily compile anything from f77 to Fortran 2008)
That would suggest that, for example, C11 should not be referred to as the C programming language, since it introduces some new features? Similarly, if I said “I have a Camaro” to you, would you expect that I have a 1970s car or something more recent?
C11 still looks like (and still is) C, despite the new features. In contrast, to say that the new Fortran is the same 'thing' as the old Fortran is disingenuous. For example, the old Fortran did not support recursion (which fact I consider a feature rather than a "bug").
ANSI C looks pretty different to K&R C, especially in the function declaration syntax. The only difference is Fortran had its facelift just a little later.
On the other hand, C was given a different name when classes were added to it.
Thing is, adding and/or removing features may have an impact on the implementation, can potentially introduce subtle changes to the semantics of the existing constructs and therefore their behavior at run time, have general performance implications (e.g. due to a change in what optimizations are possible), etc. Is it "the same language"? Not quite. It is a "super-language," at best. If it's 100% backwards compatible, that is (which is not guaranteed).
I've done a bit of HPC work, and Fortran is very much "still a thing" there. MPI[1] is pretty much the only game in town for between-node parallelism, and MPI comes with interfaces for C, C++ and Fortran. The only other language that's even been run at petascale is Julia, and as far as I can tell that's still using Julia's `ccall()` under the hood to interact with the MPI C libraries [e.g., 2].
Certainly legacy code is part of the picture, but not always as directly as one might think. Probably the biggest factor is that Fortran compilers tend to be very good -- partially a chicken and egg issue (e.g., Intel puts effort into into ifort because its HPC customers want to use Fortran), but I think there's also at some level a tradeoff between language convenience versus ease of optimizing into machine code. To give one concrete example, until C99 added the `restrict` keyword, it fundamentally wasn't possible for the compiler to optimize C code as heavily as it could optimize Fortran in certain common situations because of pointer aliasing issues.
It's probably also worth noting that modern Fortran is a long way from f77.
> MPI[1] is pretty much the only game in town for between-node parallelism, and MPI comes with interfaces for C, C++ and Fortran.
I believe the MPI C++ interface has been deleted from the MPI standard. It was a bit pointless, since it was essentially just the C interface with slightly different wording, and C++ can of course call C just fine. For people wanting a higher level C++ MPI interface, I believe the common choice is Boost MPI.
> The only other language that's even been run at petascale is Julia, and as far as I can tell that's still using Julia's `ccall()` under the hood to interact with the MPI C libraries [e.g., 2].
I don't think that's a bad thing, per se. No need to reinvent the wheel.
And, it's the same thing for Fortran really. The most common MPI libraries are implemented in C, with the Fortran binding being a fairly thin wrapper that calls the C implementation. The only real difference is that the Fortran binding is an official part of the MPI standard.
Hah, I hadn't heard that! My impression has been that good C++ written for HPC tends to look a lot like plain C anyways. And ccall is pretty efficient, so I'm not complaining there either.
That's a good point. Hierarchical parallelism is becoming increasingly important, so having one language that can be used both within-node and between-node is very convenient, and could add to the lock-in factor.
Good point and this is btw. exactly where Nvidia is heading. There will be a point in the future where you just program kernels and/or map/reduce functions and/or library functions and then call them to execute on a GPU cluster, passing in a configuration for network topology, node-level topology (how many GPUs, how are they connected) and chip-level topology (grid+block size).
The address space will be shared on the whole cluster, supported by an interconnect that’s so fast that most researcher can just stop caring about communication / data locality (see how DGX-2 works).
> The address space will be shared on the whole cluster, supported by an interconnect that’s so fast that most researcher can just stop caring about communication / data locality
There will always be people who will care because locality will always matter (thanks, physics). Improvements in technology may make it easier and cheaper to solve today's problems, but as technology improves we simply begin to tackle new, more difficult problems.
Today's chips provide more performance than whole clusters from 20 years ago and can perform yesterday's jobs on a single chip. But that doesn't mean clusters stopped being a thing.
I do think there’s a paradigm shift coming. it’s a combination of the ongoing shift away from latency- to throughput oriented design with the capabilities shown in new interlinks, especially nvlink/nvswitch. This allows DGX-2 to already cover a fair amount of what would otherwise have to be programmed for midsized clusters - if it can be made to scale one more order of magnitude (i.e. ~ 10 DGX) I think there’s not much left that wouldn’t fit there but would fit something like Titan. There’s not that much so embarassingly parallel that the communication overhead doesn’t constrain it, and if doesn’t, you again don’t care much about data locality as it becomes trivial (e.g. compute intensive map function).
If you are looking for an alternative to MPI, you should try Coarray Fortran. It supports parallel programming for both "within a node" and "between node" communication.
Coarray Fortran is now part of the Fortran programming language, and has a very simple syntax similar to array operations.
Based on my experience, the performance would depend on the compiler implementation and I would recommend GCC compiler instead of Intel compiler.
Interesting, I haven't tried that. Apparently it uses either MPI or GASNet for the between-node communication, depending on configuration? I don't know anything about the latter, but apparently it's an LBL product [1].
I used to work in HPC. The Mellanox gear, specifically InfiniBand is very good.
Fun fact: if you're working at a Saudi Arabian HPC center, say KAUST, your interconnects are purely Ethernet. Mellanox is (partially?) an Israeli company, and that's not very politically comfortable with procurement.
Better than? Not necessarily disagreeing, but I'm not sure what the alternatives even are at the same level of abstraction. I mean, there's PNNL's global arrays [1] but that's higher level, or Sandia Portals [2] which is lower/transport level.. Perhaps there are newer/alternative options I don't know about?
Global arrays is normally used over MPI anyhow. I guess there's SHMEM, but that's integrated with at least OpenMPI (and others, I think). CHARM++ has been used at scale, but it's semi-proprietary.
I happen to be one of those "real engineers", working in an aerospace company. When I have to re-do some calculations done 30 years ago (yeah, we need to have traceability of calculations for that long and more ... ) , I grab the F77 source code, build with gfortran or Intel Fortran, and believe me, it builds and runs. From Windows workstations to Linux clusters. Easy.
Then I get the python code from 6 months ago, and spend hours figuring the right python version, and library compatibilities ...
And when I read the fortran code, made by (physical stuff) engineers, it is ugly but simple - there is one or maybe two ways to do something. Not the 1000000's of ways of doing the same thing in Python. And I don't have to learn 50 new libraries to understand the code.
By the way, Python has been my language of choice for most stuff in the last 15 years ...
This was not my experience at all when trying to build Fortran code.
Which version of FFTW were they using? Undocumented. Hard binding to MKL as opposed to generic LAPACK usage? Ok, Intel fortran only then... Still missing some symbols? Oh, that's a symbol from this other esoteric library from the same group. Better see if I can download the sources for that. Etc.
There are projects without dependencies, then there are projects with dependencies in languages with sane development oriented package managers, then there are projects with dependency hell. I found fortran projects usually fell into one of the two extremes, certainly Fortran is not modern in the sense that it has no package manager which I believe most people would see as a prerequisite for a modern language.
Almost every time there will be regression test cases to verify the results,and the new builds will never reproduce perfectly the numerical results. In general is not a issue, has to be handled case by case.
A bigger issue is that modern compilers will raise many warnings, and simple inspection will find many bugs, like array or loop boundaries errors, and when you fix the code you get a different result.
Then it has to be handled case by case as well, since you have to check if that affects safety of flight.
I'm working on a Fortran book with Manning Pubs [1]. You can browse some of the projects that accompany the book here [2]. At Cloudrun [3], we're running some 900K lines of parallel Fortran in the cloud for on-demand custom weather prediction. Yes, it's very much still a thing, it's just that the application domains isn't as much in the mainstream media nowadays.
You can think of R as a wrapper around Fortran, in much the same way that Python is a wrapper around C.
Not only does building R itself require a Fortran compiler, common packages do as well. If you look inside your R package tarballs, particularly of solvers, you will see Fortran.
(Scientific computing libraries in Python also wrap Fortran, but CPython doesn't depend on it. Base R depends on Fortran.)
There is just a ton of well written Fortran for number crunching and there is zero reasons to not use them. You wouldn't gain speed and you would loss the decades of stability these Fortan scripts have given the scientific community.
I am guessing all of us old timers remember the pain of Fortran in the past.
I had to compile some Python libraries from C source on a Solaris box once. Initially I got errors due to an incompatible version of the OS maths libraries that were clearly written in Fortran, though I can't remember the details of the error messages. So not really a direct dependency on Fortran, but it was still in the mix.
In one of my classes, as a gross oversimplification, I explain that what determines the difference in speed between an R and a Python implementation of something is how much of it is secretly written in C or Fortran.
I worked with some numerical analyst/software people before. Although their C++ code looks modern and fancy, the effort to install all the dependencies and finally utilize those is not worthy of the trouble. Eventually, it's all about solving the problem efficiently, not the appearance.
I make this joke myself all the time and yet it does sting to see someone else write it.
I do generally agree with the point though, and sometimes I feel that I and people I work with become hung up on things like code elegance and testability rather than actually addressing the task at hand.
Fortran II was my first programming language, and a key element to my first summer job, and to my first FTJ. There were a lot of things about it that eventually made me dissatisfied enough that I studied compilers to see why there were such goofy restrictions, dabbling in XPL. I haven't been back, but germane to this article, I visited with my professor with whom I had the summer job 40 years later, and he pulled out the code that I had written that was still in production, It had been changed a bit, but it was still very recognizable. And it survived because, while it was a small part of a very large numerical analysis project, it worked.
The article makes me wonder if I had a problem of this sort again, would I choose Fortran. I'm thinking not, but then again, there are hints that the language has changed significantly over the last 50 years.
I have used many other languages during that time, and am fond of several of them (bliss 36, lisp, various assembler) and dislike others (RPG III, COBOL, perl). Python would be in the middle.
Tons of linear algebra and optimization code is still actively developed in Fortran. The most accurate solver I've found for small quadratic programs is a Fortran code maintained by a grumpy German math professor.
Quite a few languages people declare "dead" are doing just fine. We tend to associate languages like FORTRAN, COBOL, and RPG, with the term legacy which implies old an obsolete. Yet these are behind systems we use everyday but just are not aware of it. Plus these languages serve their purpose very well, FORTRAN is excellent in scientific solutions, COBOL and RPG excel in business math and similar.
Solutions to problems are not one language dependent. Exiting systems can have new modern customer facing solutions implemented all the while exploiting the existing code base. In fact there are times were its desirable to keep what is known to work and just give new means to access the code and data behind it.
> [D]eveloped by IBM in 1959 as the Report Program Generator - a tool to replicate punched card processing on the IBM 1401 then updated to RPG II for the IBM System/3 in the late 1960s, and since evolved into an HLL equivalent to COBOL and PL/I.
People who proselytize programming languages are missing the big picture. Its always about making computers do stuff and programming language of choice is just an implementation detail.
Having said that, we do need to look into other kinds of computers, I see none of the proselytizing there..
>People who proselytize programming languages are missing the big picture. Its always about making computers do stuff and programming language of choice is just an implementation detail.
What if you're missing the big picture?
Many developers don't care about the stuff we make computers do, but are passionate of the tools and implementation details of doing them.
If I work on some accounting enterprise software I could not care less for accounting and the associated business goals. But I could very much like working with this over that tool and finding nice ways to solve various problems towards the (still boring to me) end goal.
Unlike we're working on our personal/hobby/passion projects, the "stuff we make computers do" are usually irrelevant to us -- business owners care about them. But, as programmers, we presumably do care about programming and programming tools and concepts.
>Many developers don't care about the stuff we make computers do, but are passionate of the tools and implementation details of doing them
This is a positive characteristic of a hacker. Very much desirable.
Do note that I was specifically talking about proselytizing in a negative sense. One can love their tools and still not participate in dick measuring contest of which tool is better and actively belittle someone else who also loves their own tools.
There can be a constructive discussion about programming languages (it does happen among the PL theorists and hackers), but in general PL theorists are a rare sight.
In a general forum, sadly, many developers do judge a person with what programming language they use, disregarding what they have created. This goes against the foundation of hacker ethics.
As someone working with scientific codes in (relatively) modern Fortran, I agree that Fortran is still very strong in its niche. I don’t know any other language that makes it as easy to write fast numerical codes, and it’s relatively hard to accidentally mess up the performance.
That said, I still find it very much feels outdated as soon as you veer off its core features. Handling strings is a major pain, writing short functions feels overly verbose, and proper testing is difficult. The language is full of arcane details that require reading up on best practices even for simple things such as defining the precision of your numbers.
So while Fortran still is strong in some areas, I do see it as an outdated language, just one that doesn’t really have a good replacement yet. Hopefully Julia can replace it for most use cases, or maybe Rust if they figure out their HPC story.
> or maybe Rust if they figure out their HPC story.
Does rustc currently have the ability to implement the right optimizations? One of the big things that FORTRAN can do is prevent aliasing. Can you prevent aliasing in LLVM?
Yes, through the “noalias” attribute. Rustc uses it, well, when it’s not buggy at least. It’s currently turned off pending an upstream fix; we’ll have it back on when it’s fixed.
Fortran is still highly relevant in the science community. I worked with a couple senior scientists from a national lab and they mainly only knew Fortran and C/++. Almost all their legacy software was in Fortran so everyone who joined the group had to learn it too. Not coincidentally, they mainly teach Fortran and C++ during my physics minor (same programming courses as majors). I know they introduced Python to the curriculum these past few years though.
I find that the industry in general has a bias against gradual improvements of existing projects. Several times throughout my career, I've worked on mainline supposedly-legacy products while other teams work on greenfield rewrites, and the result is invariably that I'm able to improve the left-for-dead legacy project faster than the rewrite group can make something new --- and up with an objectively better product.
The difficulty of continuity is overstated. It's very rare that something is so far gone --- even Fortran --- that it can't be improved by insensible degrees into something just as modern and usable as something started from scratch.
Everyone I know in astrophysics or engineering who can uses something else--Python, C++, and Julia, are big winners.
The problem is that there are these huge monolithic codebases written in Fortran that are difficult to get out of using. People end up using Fortran because they need a hydrodynamics code, or a navier stokes solver, or a chemical kinetics code, but they can't justify building one from scratch--that's not how you get a Ph.D., and Ph.D. students and post docs are the people doing the brunt of this work.
There's also the added disincentive that even if you do re-implement a solver in a modern, maintainable way, then your results might disagree with previous published work--not necessarily because your results are wrong, but because the published work probably has bugs that no one has been able to discover or reason about because the code is so complex.
> Everyone I know in astrophysics or engineering who can uses something else--Python, C++, and Julia, are big winners.
I know many other languages, but for engineering work I have started some projects in Fortran simply because it is a very good fit. Python is not fast enough, C++ is a huge mess, and Julia has just been released, it is not mature enough yet.
You would be surprised of how many people are writing solvers from scratch. It is not the most common job in research, but it still is an active field.
And your point about results mismatch is true, this is a very real problem, but independent of the programming language used. To give a concrete example, many people have reimplemented almost identical material subroutines for FEM packages and obtain different results, but this happens to people using Fortran and people using C++.
What is modern Fortran good at? It seems like Fortran's niche has been scientific computing, specifically high performance linear algebra. (I wonder if that makes Fortran good for neural networks.)
In college, I remember attending a great talk by a visiting Stanford compiler prof. who talked about being able to produce faster C++ code by transpiling to Fortran first than by optimizing the C++ or IR directly. Not sure if that's still true, but at the time it seems Fortran was a more restrictive language which permitted stronger optimization.
> Not sure if that's still true, but at the time it seems Fortran was a more restrictive language which permitted stronger optimization.
There is no aliasing in Fortran. In C (and using compiler extensions in C++) you manually use restrict to give the same information to the compiler, but it is a trickier process than in Fortran. Rust inherits this benefit from Fortran, but last I checked they had all related optimizations disabled (and most numerical computing friendly features in Rust tend to be on the back burner).
That was due to a bug in LLVM when noalias was used in conjunction with LLVM’s stack unwinding facilities. It’s since been re-enabled now that the upstream bug has been fixed: https://github.com/rust-lang/rust/pull/50744
> most numerical computing friendly features in Rust tend to be on the back burner
Yes and no, we’ve been working on underlying stuff that will be needed in order to ship those features, and while a roadmap for next year has not been decided, most believe it will be a major component.
Definitely still scientific computing ("formula translator" after all). In particular it's still frequently used in the supercomputing/HPC space for large parallel numerical simulation codes.
The issue with optimizability in C is mostly solved by the C99 `restrict` keyword, but Fortran compilers are still very good for numerical work.
In practical terms: what kind of performance difference would one see between optimized C(++) and optimized Fortran when running exactly the same algorithm?
According to the Julia benchmarks (and I assume Julia folks are relatively unbiased re: C vs Fortran), it depends on the algorithm [1], but less than an order of magnitude in either direction. I'm kinda surprised Fortran doesn't do better there actually -- but at this level it probably also depends on compiler and architecture.
Thats why we still use it at my work (CFD). Its easier to produce fast matrix/linear algebra code than C++ for a scientist who is not a Software Engineer.
This has been my experience. Although I know other languages, when I go to write a solver for CFD, radiation transport, heat transfer, etc I normally have gone with OpenMPI and Fortran. For me, Fortran is almost like writing psuedocode.
Fortran produces efficient out-of-the-box code for dealing with multidimensional arrays. That said, most hand-written fortran doesn't compete with an optimized BLAS or tensor library using CPU-tuned intrinsics. For GPUs, hand-optimized libraries in CUDA rule and those are generally called from C++. Hence why you see deep learning generally in C++ (wrapped in python).
I would love to see a specific example problem, I assume it would be a numerical / linear algebra problem, that is well suited for Fortran. I want to get to the bottom
of where the speed / accuracy comes from w a Fortran compiler.
"Many believe that it’s complex and outdated. Programming in Fortran is perceived like riding a horse-driven carriage to work. Fortran?"
That's incredibly ignorant, and shows very clearly how much we suck ass as an industry, and why we're such a bad industry to work in. I'm currently working on LAPACK, and that's written in F90. Why? Speed!
To add insult to injury, Fortran is really nice to program in, especially modern Fortran, especially targetting the GPU's.
Speed and trust. Trust in the algos that were devised and implemented through the history.
Older codes often had a parameter to control a number of significant digits to maintain. Perhaps a limitation in those days, but also it allowed to control error tolerance. Indeed, engineerig it is.
Article claims that "modern Fortran already has [...] generics" but AFAIK it's still not possible to write e.g. a type-generic sort function (unless you restrict yourself to types all inheriting from a common class, a bit like very early Java). So there's no standard library of basic algorithms and data structures, and you can't even write one, without resorting to meta-programming (likely in a language other than Fortran). From what I have seen about Fortran2018, this will not improve any time soon.
You can use "class(*)" and the "select type" construct, and yes, polymorphism.
So to compare to C++ or Java, it a bit like inheritance and RTTI. Not generating type-specialized functions at compile-time like with C++ templates, which is what's needed for high performance.
Scientific computing doesn't (only) mean small tight loops doing linear algebra. Scientific applications can be millions of LOC, with all the attendant problems of managing complexity at scale. Polymorphism can well be a useful tool to have in the toolbox. Just because it's not an appropriate tool when writing high-performance kernel code, doesn't mean it's not useful elsewhere.
I quite like having functions which take generic "scalars" (double/float/std::complex<double>/boost multiprecision etc) without me having to worry about the specific type everywhere. Of course, sufficiently far down eventually one has to care about the scalar type, but on a high level, it’s nice to ignore it.
The same goes when building abstractions based on this. Tensor<T> instead of Tensor_double, Tensor_float etc.
And then you have all the typical usecases of polymorphism in the wrappers and callers and executables and whatnot which are also to some degree part of "scientific computing".
One example I often give is the nabla operator from vector calculus that is in one form or another used in any fluid dynamics code. If you define nabla as a function that operates on a Field type, and Scalar and Vector types extend the Field type, then nabla can be coded to return a gradient (Vector) if passed a Scalar, and divergence (Scalar) if passed a Vector. Polymorphism is an elegant way to implement this.
For something quantitative: the most time on general-purpose HPC systems is likely to be taken by materials science software written in Fortran. There are numbers somewhere for Archer, the UK "tier 1" system, for instance.
VASP is likely to be top. I've not seen much use of abinit. Gromacs, charmm, and siesta don't use Fortran. See also cp2k, elk, quantum espresso, nwchem, and dl_poly, amongst others, which do, and are common. Gaussian is unfortunately common, and I think is Fortran. That's on top of important Fortran codes in ocean and meteorological modelling, and engineering. https://www.archer.ac.uk/status/codes/ has a lot of "other", partly because they don't log what actually runs with whatever is in the Cray compute node kernel which would do a similar job to audit, for instance.
There are certainly much more Fortran users today then 30 years ago. It is not only widespread in science and engineering education but subsequently of course also in the research codes. However, there is a gradual shift towards scientific python when high performance is not needed. But as the author notices, Fortran is the prominent language for high performance programming because it is so simple and has linear algebra built in at its hearth.
> It stll excels in the good old structured programming. It has features that mainstream C-like languages lack. For instance, it can exit or continue a loop from a nested loop.
You can do this in C with gotos. The other things that are mentioned are similarly easy to implement–they just don't have a particularly nice syntax.
Fortran has computed goto's (though marked obsolescent). Most mainstream C compilers support computed goto's as an extension, Visual Studio being one of the oddballs.
Computed goto's can really speed up state machines by allowing you to replace conditional branching with unconditional indirect branching. I've never seen switch-based solutions come close, notwithstanding that switch statement optimizations are supposedly part of the bread-and-butter of optimizing C compilers.
Note that GCC's documentation on computed goto's uses an example where they index a static array of label pointers. AFAIU the example does this because it speeds up linking and because its more analogous to Fortran code. But it's easier and more performant to simply store and use the label pointers directly (in place of an integer for indexing the array), and the additional linking cost is irrelevant compared to the amount of symbol linking required for even a simple C++ program.
That makes sense. gotos make code generation so much easier, and computed gotos are even more powerful. Lua 5.2 added goto precisely to make it easier to machine generate Lua code. WebAssembly lacks gotos and it's a serious impediment--compiling to WebAssembly requires a special "reloop" algorithm which, while conceptually simple for most code, doesn't work well for the type of critical code for which using goto and especially computed gotos are useful for. Reloop doesn't always work and so sometimes the only way to compile to WebAssembly is to add indirection, such as compiling to an intermediate state machine.
I think WebAssembly originally required reloop because it was originally designed as a purely stack-based VM with properties that made it trivial to verify the bytecode (similar to BPF). Then things got much more complicated for performance reasons and I feel like they probably could have added goto support at that point with little marginal complexity. It's clearly possible--Ethereum did it, BPF permits forward jumps, and eBPF permits both forwards and backwards jumps.
In general, you can implement anything with any programming language. Having a particular nice syntax (and semantics) usually is what makes a language better suited than another one.
I think it has more to do with national security; i.e., not sharing our tech with other nation states that might use it against us. That used to be more of a thing a generation ago. For example, it used to be a crime to export certain computer technology to communist China. Now lot of our stuff is built there.
The reason is nationalism, what I wrote in the parentheses was meant as a tongue in cheek paraphrasing of some isolation supporter, I thought it was obvious because I would never call a person "alien" (as in intruder) seriously. The reason here is "nationalism" because that includes the whole history (which is in turn the reason for it).
FORTRAN is definitely still used all over the place in scientific computing... If you go to a CS department that does some engineering work and work on numerical you might have to write some.
Most languages can access C without much hassle (Ruby, Python, ...), F77 predates this, and while you can combine C and F77 for a particular compiler, it is a painful and unrewarding task (and the result is compiler-specific). It is often quicker to run the F77 though f2c so you have an entirely C codebase (even if f2c is not really readable).
Both BLAS and LAPACK are written in fortran aren't they? So that's all the core numerical routines for numpy/scipy and R and I assume many other things.
Sort of. The bits of BLAS that make it performant aren't written in Fortran -- either assembler, or the equivalent SIMD intrinsics. (See OpenBLAS and BLIS source.)
It recently popped up in a place I didn't expect. In checking out some of the deep learning work around neural style transfer, I found that many of the examples make reference the L-BFGS optimizer. Peeling a layer back, it appears that a commonly used implementation of L-BFGS is in scipy, which is a wrapper around a Fortran implementation.
if Fortran is built for scientific computing, how does it compare to modern languages that are used for scientific computing for the same area like Julia, R, Python?
What is one of the first statements made in every software license ever written in the history of computing?
Civil Engineers: Bridge falls down, PE's stamp is on designs, design faulty, is liable = Real engineer.
Software "Engineer": Program fails, money lost. "LOL must be bugz, btw: THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM “AS IS” WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION." = not an engineer.
Can you imagine if an aerospace engineer submitted the designs for a new turbine engine with an asterix at the bottom of the blueprints saying "this engine design is provided as-is without warranty of any kind either expressed or implied, including, but not limited to whether or not it will blow up in flight"?
The programmers for the Space Shuttle program may have been actual engineers, because I'm willing to bet that they submitted documents guaranteeing the performance and operation of the code they wrote.
Yep:
>NASA knows how good the software has to be. Before every flight, Ted Keller, the senior technical manager of the on-board shuttle group, flies to Florida where he signs a document certifying that the software will not endanger the shuttle. If Keller can’t go, a formal line of succession dictates who can sign in his place.
So maybe some programmers can be engineers, if they put their money where their mouths are and don't go "lol bugz, didn't buy support contract" if something goes wrong.
It'd be great if "software engineer" == "real engineer" (well, for some - for me it'd suck, because I don't even have a good degree, let alone any kind of background or education for engineering certs) - but imagine if that were true...
...you could basically kiss every kind of cheap software, home computing, mobile phones, etc - goodbye.
Because nobody would be able to afford any of that. No one in business would pay the money for thorough design and vetting of business software, let alone consumer-grade code. Ordinary consumers would find any such software prohibitively expensive (on the order or worse than trying to buy a license for a commercial Unix back in the 1980s - assuming you even had the machine to run it on).
Heck, the amount of time alone it would take to carefully design such systems, then build, debug, test, verify, etc - it would be insanely and prohibitively expensive.
I'm sure you know this.
Which is why you only ever see this kind of effort applied to software used in aerospace and medical fields (and probably to a lesser extent automotive, as well as commercial software used for engineering - like CAD/CAM, FEA, etc).
Even there, costly (sometimes deadly) mistakes can and do happen.
Regardless, though - the cost to create that software is phenomenal (and the cost to maintain it even more so - every single minor or major change has to be carefully considered and researched, then if put into place, carefully tested and vetted to ensure nothing else is broken by it).
It ultimately comes down to money and the market; if all software were held to the standards of say aerospace software, nobody could afford it, nobody would sell it, and nobody would buy it. The fact that only a few particular industries have settled on such reliability being important says it all. Had all industries and concerns wanted such reliable software, we would see it today. Since we don't, it is likely because it was deemed too expensive for little gain overall.
You can certainly blame the builder, but the "real" engineers I've talked to all spoke of traditions like a _ring on their pinky_ so that when they are signing their name to a design, they are reminded (as the ring touches paper) that their signature affects real people's lives, and that both their professional reputation and personal liability might be on the line if they screw it up.
If an engineer designs my bridge, and it collapses when more than one car is on it, he can't laugh at me for not giving him a spec that included that requirement. Contrast this with software, where if a database stores plaintext passwords, the programmers are probably protected if there is sufficient documentation that the bosses/customer insisted on it.
I'm grateful to have the flexibility of creating effectively magical things that can change easily to meet changing requirements, but I also feel like there will come a time when we will need to have some similar degree of engineering rigor mandated for software.
It's beyond traditions in some regions as well. For example, in Canada, the "engineer" title is protected, so anyone using that title must be a registered member of the provincial association (wherever they practice).
Do you think an engineer just sits down, draws your bridge up after some math and then builds it? It takes a team of civil engineers multiple iterations with tests, setbacks, feedback loops etc to get to the final design.
That's actually quite rare and certainly in some engineering disciplines no way would you wear a ring whilst working for H&S reasons I am thinking of EE and Mech here.
Also professional engineers (who are the only engineers able to use the title "engineer") do P.Eng exams, and must maintain education / practice credits throughout their career. They must also practice as an Engineer in Training and record several years of experience before even applying. So yes, a "real" engineer is quite a lot more involved than new CS grad getting with 1 year experience getting a job as "Senior Software Engineer".
Years ago, I had an amusing fight with the C++ committee over this.
In C++, you can overload the "[]" operator:
But you can't do that for multiple arguments. Why? An obscure feature of C and C++ is that the precedence of the comma is different for function argument lists and subscript lists. In is a function call with two arguments. is an invocation of the comma operator, which ignores the second operand. The syntax charts reflect this. So I proposed that the syntax for comma usage be the same in both "()" and "[]" forms. If you see "foo[a,b]" in a C/C++ program, it's almost certainly an error. I searched through large volumes of code and was unable to find a single use of the comma operator inside brackets.But no, there's a "save the comma operator" faction.
C people just do not get multidimensional arrays.