I keep looking around and not finding any, so let me just try here before someone just takes it and slopifies it:
* built-in multidimensional arrays with efficient storage.
* related to this: built-in array intrinsics
```
real, dimension(100,100) :: A, B, C
C = A + B
```
this kind of code is already a close-to optimal "naive" implementation (not considering parallelization). so you start already at a solid place. then you can easily run it in parallel without too much specialized knowledge with OpenMP, OpenACC, MPI or even CUDA. the only thing you really need to be aware of when implementing your own loops/kernels: the intrinsic storage order, to optimize for cache hits.
* crucial: all the above amounts to a standard/best practice about how data is structured and formatted. everyone just uses the built-ins. Thus, interoperability between native Fortran numerical libraries is usually a complete non-issue. Meanwhile, Cpp has a fractured ecosystem with different array/vector types for its libraries. Converting between one and the other is usually a no-go.
* next, the intent plus pass-by-reference system. it combines IMO the best of both worlds of a functional vs. procedural approach:
- I can predefine if an array / variable is intended as input, output or both, simply with `intent(in)`, `intent(out)` etc.
- compiler can thus do checks for me if I'm breaking a contract.
- yet, since I pass by ref, I don't have to worry about memory not being used efficiently. This really matters once you deal with GB, TB or even PB worth of data going into a simulation over time - there's just no way to deal with that in a purely functional way.
- only where really necessary, you can still drop down to pointer semantics, e.g. for the outer glue code of a simulation that swaps outputs and inputs for the next time step.
- having `restrict`-by-default semantics here helps immensely. Imagine you have many arrays with lots of data to deal with, and your kernels access always several at a time to do calculations. In Fortran I can write it intuitively, while in C/C++ I must remember to specify `restrict` or to preload all point-inputs first into separate variables to direct the compiler. Otherwise the memory pressure increases, and by far most such physics simulations are memory bandwidth bound - every access counts.
* finally, a clean symbol definition system that decouples types from byte lengths. a `float` in fortran is just `real(4)`, a double is `real(8)`, a long int is `integer(8)` and so on. now, it's trivial to do a bit of preprocessing to switch the precision.
However, the last part is where Cpp has a strong advantage: Well supported meta-programming (generics, templating or even just well supported pre-processors). Fortran's compilers come with a lot of built-ins, so the lack of these is less of an issue than you might think, but it's still a limiting factor. All that being said, a typical scientist doesn't tend to care and just wants to solve a particular problem rather than thinking in generalized frameworks - and that's why I find Fortran still serves them better for numerics than anything that came since.