Is Fortran easier to optimize than C for heavy calculations?
stackoverflow.com
stackoverflow.com
For most of the class, FORTRAN felt mostly like a penalty box. It was full of behaviors and limitations that dated straight back to the early 1960's, if not earlier. It was at least easy to make it do the things it could do.
Where my opinion of FORTRAN changed was in the last time that class met. By that point, we 'knew' the language and the instructor was on to more advanced topics, including optimization. He made a point of how an optimizing FORTRAN compiler could re-nest three nested loops for the purpose of automatically improving memory access patterns. Vector Cray machines (which implemented common looping patterns in hardware) shipped with FORTRAN compilers that could take normal-looking loops and turn them into code that optimized for the hardware.
I don't know if it's still the case (thanks to restrict and better compilers), but he made a compelling point at the time for the simplicity and restrictions of the language leading to more options for optimization.
Nobody (really) uses restrict, certainly not at scale.
This is especially obvious because when Rust came along (which can effectively use restrict in most places, much like Fortran) it uncovered a lot of serious bugs in LLVM's implementation for restrict. That whole situation only really stabilized in the past year or so.
Maybe GCC does a better job as a result of having gfortran.
It is the primary language for some of the most intensive super-computing tasks, such as in astronomy, climate modeling, computational chemistry, computational economics, computational fluid dynamics, computational physics.... Since the early 2000s, many of the widely used support libraries have also been implemented in C and more recently, in C++... For this reason, facilities for inter-operation with C [1] were added to Fortran 2003 and enhanced by the ISO/IEC technical specification 29113, which was incorporated into Fortran 2018... " [0]
[0] https://en.wikipedia.org/wiki/Fortran#Science_and_engineerin...
[1] https://en.wikipedia.org/wiki/Foreign_function_interface
* Fortran arrays do not alias by default. (This is the most well-known, and is indeed repeated in the answers on StackOverflow)
* Fortran N-D array indexes are UB if they go out-of-bounds in a way that isn't true for most other languages. In C, if you have a declaration `int x[2][3];`, the expression `x[0][4]` is a legal way to access one of the elements--it's only UB to go outside the bounds of x as a whole, not individual rows. In Fortran, the equivalent expression (well, modified because Fortran is column-major) is illegal.
* Fortran has built in array expressions--you can express a vector addition just by adding the arrays themselves. This makes it easier to autovectorize.
Is it really? Does the C-compiler guarantee that your 6 ints will be stored in a contiguous memory area?
int *x; int numRows = 2; int numCols = 3;
x = (int *)malloc(numRows * sizeof(int ));
for (int i = 0; i < numRows; i++) { x[i] = (int )malloc(numCols * sizeof(int)); }
// This could segfault x[0][4] = 42;
Even weirder, for that declaration
int x[2][3];
I think the following are legal expressions: 0[1[x]]
1[x][0]
Or is there something in the language that makes that break parsing?Even though multidimensional arrays are laid out contiguously in C/C++, it doesn't follow that indexing can traverse from one row to the next. By C17 6.5.6/8, pointer-integer addition can only create pointers to other elements (or one past the last element) of the same array object that the original pointer came from. However, an int[2][3] array has two int[3] arrays as elements, and each of those arrays has three int objects as elements. The int objects are not all elements of the same array (the standard's definition of an array element, 6.2.5/20, doesn't apply recursively), so you aren't allowed to index between them with expressions like (x[0] + 4).
Compilers seem to be pretty lenient about this in practice, though, only optimizing out out-of-bounds accesses when they exit the complete object. But at least Clang's UBSan will complain if you try to index between rows [0].
A Fortran program that reads its input, calculates and finally writes out its output does not have to execute any particular instruction at all. As long as the answer is "AS IF" it had done the user-specified computation, the Fortran compiler has done its job. In between I/O, it submerges into the ineffable like a Cold War SSBN.
C is about the instruments, the players and the conductor, Fortran is about the music.
See "C is not a low-level language" - https://queue.acm.org/detail.cfm?id=3212479
C, as an operating system implementation language, is trying to do something fundamentally different than Fortran.
You live by memory address, you die by memory address.
This would also be true for assembly, hardly a high level language
I expect the hardware to handle cache coherency in that situation. What the compiler does should be irrelevant.
So maybe it's accurate to say that C is "more compatible" with real hardware, in the sense that its abstract machine is more isomorphic to what's really happening than Fortran's is. But it's not exactly "closer to hardware" in the way we might be tempted to think; it's more of a lingua franca that your processor happens to speak.
If you're still tempted to consider C "close to hardware", consider that you can compile the same code for a Z80 and a Threadripper. What hardware exactly are you controlling that's common to both?
As PhilipRoman said, this is also true of assembly (or any other programming language model[1]).
> If you're still tempted to consider C "close to hardware", consider that you can compile the same code for a Z80 and a Threadripper. What hardware exactly are you controlling that's common to both?
In both of them I can write to a memory-mapped I/O device, if it has one. I can write a custom memory allocator for a pool that I'm managing myself. I can't do either of those in Fortran or Javascript.
[1] Why does it have to be true of any other programming language model? Well, maybe I exaggerate slightly. But can you show me a (single threaded) programming language where "a = 1" does not mean that on the next line, a will be 1?
Generally agree with your point, but just to play the devil's advocate, in a CPU with exposed pipeline and no interlocks, setting a register to a value doesn't guarantee that a following instruction reading from that register will see the last value written.
MIPS I.
https://news.ycombinator.com/item?id=30022022 How ISO C became unusable for operating systems development
Fun facts, Walter the original author of D language wrote his popular Empire game in Fortran [2]. Some of the ideas that make Fortran fast is incorporated into D language design and this makes D is as easy if not easier to optimize than Fortran [1].
[1]Numeric age for D: Mir GLAS is faster than OpenBLAS and Eigen:
http://blog.mir.dlang.io/glas/benchmark/openblas/2016/09/23/...
[2]A Talk With Computer Gaming Pioneer Walter Bright About Empire:
https://madned.substack.com/p/a-talk-with-computer-gaming-pi...
https://gist.github.com/rbitr/3b86154f78a0f0832e8bd171615236...
if (dims in range A)
gemm_implA
else if (dims in range B)
gemm_implB
else if (dims in range C)
gemm_implC
// several more implementations...
The immediate call function is invariably a switch over several implementations, and by the time you reach that call, you have actual knowledge of dimensions.There is no universal optimal algorithm.
But yes I agree, at some point it's always going to be back to something more customized.
my professor who teaches computer programming (who previously worked as a nuclear physicist) loved mentioning Fortran whenever possible.
He truly enjoy the language and believes it is well suited for heavy computation. He said its used a lot for weather forecasting.
He is quite old now so I don’t know how true it is today, but he made me more interested in Fortran and I might give it a shot instead of C sometime soon.
It is much, much better than C for all these applications. Nobody in the field seriously considers C; all the non-Fortran codes are in C++. Both the Fortran and C++ bits have Python wrappers anyway.
Anyway the urban myth in the IT dept was, they had a custom fortan compiler written in early 70s that saved the firm from buying 2 additional IBM mainframes with the compiler writer compensated well
Fortran >= 90 has ALLOCATABLE variables (and now components), which are dynamically managed objects that are akin to std::unique_ptr<> or std::vector<> smart objects. A compiler can assume that an ALLOCATABLE object is free from aliasing, unlike a POINTER.
Fortran does have POINTERs, but they can only point to objects that have been explicitly tagged with the TARGET attribute. Objects that are not TARGETs are safe from aliasing with POINTERs. Fortran does still not have the concept of a pointer to immutable data.
Do any of these still yield performance advantages over C/C++? Probably not.
This is one of the legacies modern C++ is still wrestling with today, supposedly solved by Rust, Zig, Val and the likes.
c is call by value, you have to at least try to alias
While it may be true that the rules aren't generally checked by compilers, Toolpack did it for pure Fortran77 in the 1980s, though I don't know how comprehensive it actually was.
The technicality is that Fortran is call by descriptor (not by reference, not by value). That's a powerful difference but a sharp two-edged sword.
Sometimes (most of the time, to be fair), the compiler can work out the best solution. Sometimes, the developer has to jump through some hoops to avoid unnecessary copies.
This is only very slightly exaggerated.
In practice, it took years for Rust access these optimizations because the noalias annotations kept hitting codegen bugs in LLVM. It was an area of the compiler that had not seen much active use before Rust, and they were causing a lot of dormant bugs to become issues.
The feedback cycle was also very long because LLVM didn't fix the problems quickly. I didn't pay close attention so I don't know if this is prioritization or because those fixes required major surgery, but Rust sometimes had to wait for major version bumps in LLVM before seeing if "this time was the last time."
How much the optimizer can take advantage of that is another matter, due to what you said. Doing a quick translation to more idiomatic Rust[2], it does hoist the accesses to `matrix` out of the hot loop, and does also seem to be a bit less moving stuff around compared to the C version, which I think is putting the `matrix` values in the right places, which Rust did before the loop.
[1] fn transform(output: &mut [f32], input: &[f32], matrix: &[f32], n: &i32);
Normal C/C++ developers doing normal C/C++ development in GCC probably never saw any of those bugs.
I haven't kept up with the current state of things, but if you can finally liberally use noalias in LLVM without having to fear miscompilations, you probably have solely the increased significance of Rust to thank for that.
The #pragma GCC ivdep annotation also enables this in loops. C23 also has the [[unsequenced]] attribute, which supports pointer arguments for optimizations like these, unlike [[gnu::const]].
Pretty much any vector|matrix operations will see row against column (or colum against row) operations in equal proportion - what you gain on the swings you lose on the roundabouts whichever way you choose to orientate storage.
Indeed, given large sparse arrays, language specific default array storage methods might not even be used - and here C is (arguably) more flexible at building task specific data structures - a column(?) of pointers to row fragments (perhaps).
The big advantage early FORTRAN had was baked in anti alias assumptions and a truckload of pre existing relatively accessible NIST | NASA grade quality numerical algorithms that had been crawled over, eyeballed, tested, and pounded in simulations by many qualified applied scientists.
Yeah. That Fortran library that has a solver for sparse complex matrices that won't produce nonsense if you give it a stiff problem, that has 50 years of people pounding on it to make sure that it's rock solid? I'm not sure that there's a C/C++ equivalent for it. I absolutely cannot write an equivalent. Yes, I can write a matrix solver that takes complex data. No, I can't write one that is as efficient, as numerically stable, or as bug-free as what already exists in Fortran.
They keep the fingerpokkens off the blinkenlights.
It's a whole other ballgame when the problem domain demands near raw pipelined handling of big data streams that don't align with any type of default internal Fortran data storage.
The classic example would be late 1980s | 1990s processing of satellite image data from ground stations - the pipelines were fat (relative to the hardware of the day) and the data arrived in BIL (Band interleaved) and|or not BIL, and all manner of jumbled formats that made sense on the aquisition side.
The tools written then, back in the day, had to gulp in big end | little end data, 8 bit, 12 bit, 24 bit, and 32 bit data, and to performantly reproject and geomontage on the fly ended up with the same forms of optimizations duplicated across a crazy circus of incoming raw storage combinations.
Whether tackled via C or FORTRAN just simply reading in data and hoping it matched native FORTRAN array structures wasn't an option.