(I should preface this by saying that my claims are all based on what I see in a specific scientific community and US laboratory environment. I can't say for sure my observations are true more generally.)

I think what you say is true in principle, basically on all counts, but Fortran's niche that its advantages serve best has been continually eaten away at by the Python/C++ combination.

The "I don't know pointers" crowd steers toward Python, and the "I care about CPU cycles" crowd steers toward C++, and that relationship is symbiotic.

Julia promised to be a general-purpose language, both as fast (or faster) than Fortran/C++ and as easy (or easier) than Fortran/Python. But it doesn't seem to have panned out as they'd perhaps have hoped.

Looking at Julia Computing industry partners, it isn't that bad.

Fortran still has a big fanbase and people actively fighting to keep it alive and thriving as much as possible. Maybe the crown will remain with C++ (in terms of fast languages), and maybe Rust will eat some of its lunch, but there seems to be a resurgence in the Fortran community to make sure it remains an important language in this space.

Julia has something similar but to a lesser degree. Its main problem has been the developer experience - the latencies (slowly getting better), lack of a great editor/LSP experience (the VS Code plugin is decent, but not yet great), lack of accessible how-to resources (people are starting to write more of these), and things like that.

I don't think Julia has been around quite long enough to make that statement yet.

I've just noticed it's now available on the system I work with, which puts it in the frame for a possible (probably far out still)future move for a massive C/Fortran codebase I maintain.

My general default if I'm going to update something is C++, but from what I see of Julia it's at least worth a look.

How has it “not panned out”? Julia has been in a useable state only for a few years, and is still rapidly developing. Prior to Julia 1.0, I would have said Fortran (which I have 15 years experience with) was the most productive language for scientific computing. Today, Julia is, by a wide margin. In fact, our benchmarks show Julia matching or beating Fortran’s performance, to say nothing of the mich increased flexibility.

In my experience "very difficult mathematical code" stems from employing little or no abstraction of the problem domain, mashing together numerical algorithms with physics formula. Using Fortran is IMO as much a prerequisite as a result of this.

I fully agree with you that scripting languages are a great fit for scientific computing. But I just don't see how this is a case for Fortran.

R or Python are still much easier to use and give you the same access to high-performance numerical accelerators via BLAS/LAPACK backends. And if you need a low-level language (e.g when developing an optimised routine that solves a specific problem), I would use a high-performance widely supported systems programming language like C++ or Rust.

Fortran seems to be around mostly for legacy reasons. And I am not sure that legacy considerations alone make a good foundation.

BLAS/LAPACK are written in fortran was the point I think. There are eigen and nalgebra certainly, but neither exports as cleanly to python, and neither is significantly (if at all) faster.

CuBLAS on the other hand makes a lot of sense to use from python if you have the right hardware.

The reference implementation of BLAS that you should absolutely not use is written in Fortran. The actual implementations are written in C/C++/ASM.

That's not to say you should use C yourself, memory management is a lot nicer in Fortran, for multidimensional arrays. You can stack allocate arrays inside functions and modern fortran has move semantics. Arrays are contiguous even if they have multiple index. It's generally a nice language for math, has good C API integration (both as caller and callee) and compiles down to small libraries without runtimes.

Thanks! I did not most any of that. Except for that it is good at math.

Loops in Python are excruciatingly slow, for no reason. Some algorithms are naturally written using loops, and Python makes this very uncomfortable, forcing you to rewrite your algorithm in a non-natural way. I hate writing numerical code in Python due to that; it just feels uncomfortable.

For scientific computation I just need multi-dimensional arrays of floats and fast loops. Fortran has both of these things natively. Python has neither.

I prefer vector types for this kind of stuff, as they are more compact than loops and naturally lend themselves to auto-vectorisation. In the end it comes down to the programming style preference.

A lot of the code used in scientific computing cannot be vectorized because the latter loop is dependent on the former loop (like time stepping in a simulation). It’s not really a programming style preference.

C++ and Rust don’t really have good semantics for scientific computing (linear algebra). They’re great for many other domains like “system’s programming”, where Fortran doesn’t compete, of course. The only language that’s similarly (in fact, more) expressive as Fortran for scientific computing today, as well as holding up for performance, is Julia

> the "winning" strategy is to have high-level scripting languages where you can ignore most of the details

From personal experience, I don't believe this works in practice. Invariably, some level of logic is implemented in said scripting language, and it becomes a huge bottleneck. Counterintuitively, when you have a central thread dolling out work to many subroutines, the performance of that thread becomes more critical as opposed to less.

"Invariably, some level of logic is implemented in said scripting language, and it becomes a huge bottleneck. "

Not necessarily and if it does, then the solution is to add that logic into the low level libary. So yes, you need at least one guy who can do this.

At that point, you’re not really using the high level language.

If you aren't using the scripting language enough for it to get in the way, you're probably not using it enough for it to make things easier either.

Just to be clear, I'm not saying it's impossible to strike the balance where it makes sense. I'm just saying I've yet to see it work, and I would, therefore, not assume such a balance is given by default, or even easy to achieve.

I agree with you. There’s also the additional cost of maintaining more than one language which makes the balance that much harder to do. I could see utility in the high level language to do something with the result if it’s not computationally intensive.

The theory is that 99% of the code is in the high level language, and 1% is in the low level language.

In my experience, that means that 80-90% of the developer effort goes into getting the FFI bindings to not be terrible. (In fairness to things like python, most of my FFI experience is with Java, which is notoriously bad at such things)

Have a look how games work. They are very performance critical - you have to calculate your simulation and draw it 60 times per second. That only works, if you have high performance low level libs. It is very common to mix C++ with lua. So yes, it is about balance, but it is standard procedure by now.

And with wasm (or even webGPU) one can do now the same with javascript and the web. Some logic in scripting is fine, but the heavy calculations should be lower level if one is aiming for performance.

    From personal experience, I don't believe this works in practice.

From personal experience, this works in practice.

It works well for web or other server-side tasks where your database is doing the "heavy lifting" computationally.

I believe you, though, when you say that you've seen failures of this paradigm.

It's pretty common to see issues where an inexperienced programmer is doing a bunch of data manipulation in 200 lines of Python/Ruby/PHP/whatever that could have been accomplished more easily, more performantly, and with more atomicity in 10 lines of SQL. I wouldn't say that means it's a busted paradigm. I don't consider bad implementations of a particular paradigm to equal refutations of that paradigm.

I do appreciate that there are probably arenas (scientific computing?) where maybe this approach doesn't work at all.

> the "winning" strategy is to have high-level scripting languages where you can ignore most of the details, which call into hyper-optimized, high performance libraries. For instance, when you're using Scipy, you call into Fortran and C almost interchangeably.

That gets you about 70% of the way there, but the remaining 30% is tricky. Unfortunately, the library-call stage doesn't easily (with generality) let you fuse operations. For example, consider the example of a finite-volume code, where the flux is some nonlinear function of the interface value, taken to be the average of left and right cells:

  flux = nonlinear_flux_op(f[:-1] + f[1:]) # Compute interface fluxes
  f[1:-1] += dt/dx*(flux[:-1] - flux[1:]) # Update cell volumes based on flux differences
  # Do something special for the boundaries

Each of the substeps here is fast. Sub-indexing these 1D arrays is trivial, and their local operations are highly vectorized. The nonlinear flux operator operates element-wise, so it's easy to compile it into a version that properly (CPU) vectorizes (and/or it might already be provided, like np.exp or *2). This code is also very readable.

However, real performance would require fusing operations. There's no need to compute a full temporary array for the interface values, but that's exactly what would happen in naive NumPy code. Likewise, there's no need to compute all of the fluxes before taking their differences.

A hand-optimized code would perform the full computation in a single scan of the array, keeping all of the data in CPU cache (blocking by a cache line or so in order to use vector registers effectively). Even if given the code naively written, an optimizing compiler would try loop fusion and get most of the way to the hand-optimized result.

Code such as this is relatively FLOPS-light compared to the amount of memory accessed, so effective use of cache can lead to an order of magnitude speed-up.

With a fair amount of work you can still achieve this result in Python, but to do so you need to use one of the optimizing Python compilation tools like numba. Unfortunately, this is not easy for the not-a-programmer domain scientist, since working with numba means understanding the intersection of Python and numba's imposed compilation rules.

ChatGPT, please write a song to the approx tune of "All you need is love" about numpy.linalg.einsum

---

(Verse 1) Gather 'round, fellow coders, lend an ear, A tale of performance, crystal clear, It's about `numpy.linalg.einsum`, you see, A bridge to speed that sets your code free.

(Chorus) Numpy linalg einsum, unlocking the door, To performance realms we've hungered for, With C and Fortran, it dances in tune, Efficiency's treasure, like the brightest of moons.

(Verse 2) Python loops and arrays, they have their grace, But when it comes to speed, they're in a different space, Einsum brings us closer, to libraries so bold, C and Fortran's might, where performance unfolds.

(Chorus) Numpy linalg einsum, a liaison so grand, Linking Python's elegance with performance's land, From indices to memory, it bridges the span, Oh, the power it wields in a programmer's hand.

(Bridge) Clever indexing and Pythonic style, Are wonderful tools that make us smile, But beneath them, a world of raw speed resides, Numpy linalg einsum, where performance abides.

(Verse 3) Linear algebra's puzzles, intricate and vast, With einsum's prowess, we reach new heights so fast, Python's beauty combined with Fortran's swift glide, Numpy linalg einsum, where two worlds coincide.

(Chorus) Numpy linalg einsum, the key to the gate, Where performance and elegance conjugate, Let's embrace the speed that we've dreamed of so long, With C and Fortran, together we're strong.

(Outro) So heed this ode to efficiency's call, Numpy linalg einsum, the ace of them all, In the realm of coding, where performance takes flight, We cherish the union that powers our might.

Numba is quite restrictive... to someone comfortable in C++, coercing my code into a form Numba likes sometimes feels no easier than writing an extension using Eigen types and pybind11.

> "the "winning" strategy is to have high-level scripting languages where you can ignore most of the details, which call into hyper-optimized, high performance libraries. For instance, when you're using Scipy, you call into Fortran and C almost interchangeably"

That's nice as long as it works, but in my experience gets ugly when it doesn't. While we can currently witness this being a successful model, the impedance mismatch between both worlds is often a pain.

I bet on a language like Rust as the winning strategy in the long run - high level like Python but still close to the metal and all in one language, one compiler, one ecosystem.

Rust shares a lot of performance advantages around pointer aliasing with fortran. In fortran, pointers can't alias by default, and in rust, the borrow checker should be able to infer that a shared and exclusive reference don't alias, allowing the compiler to assume that it's not writing to input values of a loop, etc.

In C, pointers can alias by default. The last time I checked, that made these sorts of optimizations harder. However, I'm not sure how measurable the performance impact is these days.

Fortran pointers must be assumed to alias. It's Fortran's dummy arguments that must not alias.

That kind of sounds like Mojo.

> the "winning" strategy is to have high-level scripting languages where you can ignore most of the details, which call into hyper-optimized, high performance libraries. For instance, when you're using Scipy, you call into Fortran and C almost interchangeably.

Well, no. This is python's strategy. Doesn't make it the winning strategy. Python implicitly forces multiple languages upon you. A scripting one, and a performance one. Meanwhile languages such as Julia, Rust, etc. allow you to do the work in a single (fast/compiled) language. Much lower cognitive load, especially if you have multiple cores/machines to run on.

Another point I've been making for 30+ years in HPC, is that data motion is hard. Not simply between machines, but between process spaces. Take large slightly complex data structures in a fast compiled language, and move them back and forth to a scripting front end. This is hard as each language has their own specific memory layout for the structures, and impedance matching between them means you have to make trade-offs. These trade-offs often result in surprising issues as you scale up data structure size. Which is one of the reasons that only the simplest of structures (vectors/arrays) are implemented in a cross language scenario.

Moreover, these cross language boundaries implicitly prevent deeper optimization. Which leads to development of rather different scenarios for code development, including orthogonal not-quite-python based things (Triton, numba, etc.).

Fortran is a great language, and as one of the comments pointed out, its really not that hard to learn/use. The rumors of its demise are greatly exaggerated. And I note with some amusement, that they've been going on since I've been in graduate school some 30-35 years ago. Yet people keep using it.