I am getting a bit confused in this thread about what exactly is the current discomfort? Because amazingly there are multiple triggers went off. Facetious but nevertheless intriguing. Is it me saying the past is incompetent so it tickled certain folks or is it the fact that old code having bugs or me being passive aggressively accused of lying/exaggerating about the number of issues?
Let’s first clarify; maybe my english is not good enough but this is what incompetent means. Incompetent Definition & Meaning - Merriam-Webster Pick any one of them as what I meant. So I don’t care about your feelings about the past or your reminiscence of the olden days. Fortran code from 90s, (just like any code from the 90s in any language) is not capable of doing what current day code is expected to do. So let’s not insult each other’s intelligence. Users on this forum, openly ignore this fact however the code has to run on many architectures without causing trouble. Look at our CI/CD pipelines how many architectures we have to release SciPy for. So the language and the codebase should have been modernized but it didn’t. Not my fault. Not yours either.
So I hope noone here is going to attempt to say, gfortran
+ some flags is up to production grade. Then the blame is on you, not me, for being simply out of touch. Plus, just because you code in fortran does not nearly imply that your code is more scientific or you have the high ground. Hence I am not going to take anyone serious if they are writing double loop for matrix multiplication. Let’s agree to disagree on that. It doesn’t matter what anecdote you might have. That ship has sailed. And I am not going to reply to anyone on this. If the “scientific” language doesn’t give matmul
infix operator out of the box then it’s not my feelings that need to get hurt. Example from internet fortran - Why does a manually programmed matrix multiplication combined with matrix addition give better performance than the intrinsic functions? - Stack Overflow That is an example of incompetence. I am not going to accept writing 20 lines of code to do (A * b)' * (C * X)
in matlab notation. In Python it is much uglier (A @ b).T @ (C @ X)
but still 10ish characters but unfortunately we can’t influence the mother ship that defines CPython from SciPy. However if matmul(A, B)
segfaults in a language and they don’t fix it for decades; that IS incompetence. Just to disclose, I wrote my PhD thesis on Linear Matrix Inequalities for large control theoretical problems. So excuse me but I know what a long piece of linear algebra code means. I have all the scar tissue so save your loop arguments from 70s to the uninitiated.
Anyways, since the question was about bug reports, let’s take PROPACK. This piece of code has served millions of problems to date. It is, not only being my favorite old but gold code, but also, compared to other F77 code lying around, the most careful one to adhere “modern” sensitivities such as overflow/underflow detection, an effective split of concerns and a clear separation of computations, architecture awareness, OpenMP support, allowing for drop in replacement mechanism and so on. I really do think, it deserves all the credit it accrued over decades. The author clearly is not indifferent to his shortcomings of the day and left a lot of facilities for maintainability. If I’m not mistaken, that’s what he continued to do by working for Google later on.
However! this does not make it immune to the erosion of time. Let’s start with the amount of work it takes to get a decent library into a Python package. I am always smiling when I see some members of this forum saying how easy to wrap almost trivial effort to include a fortran lib because slap some iso_c_binding
here and there and off to races, well let’s investigate shall we? Take a look at the files tab of the PR.
So then after we did that, here are randomly selected things we fixed over time; emphasis on random ! so don’t get the funny idea that these are all there are. You would be very wrong
BUG: Seg. fault in scipy.sparse.linalg tests in PROPACK · Issue #15108 · scipy/scipy · GitHub (after this we made PROPACK opt-in because segfaults didn’t stop, still open and we don’t know how to fix it yet)
BLD: more PROPACK fixes, removing timer code by rgommers · Pull Request #18421 · scipy/scipy · GitHub (second
statement of fortran is both inaccurate and also anachronistic, fortran is bad in every kind of professional logging scheme, if you disagree then I need hardcore evidence compared to python logging library)
BLD: fix two `-Duse-g77-abi` regressions and a PROPACK bug by rgommers · Pull Request #18426 · scipy/scipy · GitHub
BUG: some tweaks to PROPACK f2py wrapper and build flags by rgommers · Pull Request #18263 · scipy/scipy · GitHub (read the comments in the code)
ENH: add wrappers for [zc]dotc functions for MKL compat by mckib2 · Pull Request #3 · scipy/PROPACK · GitHub (we have to fix czdotc shims because fortran/blas)
And these are just the fun ones.
In case you have been paying attention, this is a single line in my Fortran enumeration META: FORTRAN Code inventory · Issue #18566 · scipy/scipy · GitHub.
So what did I do? I read the friggin’ code and rewrote it. This should not have been me doing this. This should have been one of those lovers of fortran and linear algebra folks. But no.
There are many other lines that are causing issues or we cannot answer the feature requests; here is a recent list for slsqp META: `slsqp` optimizer · Issue #19130 · scipy/scipy · GitHub.
So what did I do? I read the friggin’ code and rewrote it. This should not have been me doing this. This should have been one of those lovers of fortran and optimization folks. But no.
Almost every line in my list has 10s of issues that we fixed here, removed there, and morphed etc. So when I say neverending I mean it.
The point here is to understand that we are not able to close these ones. That’s the problematic part. We are not a service portal number of issues don’t measure anything. We only can advise folks that we are aware of the problems and that’s about it. We don’t know how to fix in certain cases. If you ship a package for millions, there are lots of cases that your ifort and left mouse button convenience cannot provide an answer.
So another issue, compiling this f77 code for others; back in the day someone from here kindly offered to help us with NAG compiler for example over scipy. Then NAG declined to help, ifort is not free, so we have a single option that is gfortran
. I’m trying to be positive here so let’s not talk about gfortran
.
OK fine Lfortran
, I get it, it is growing and getting taller by day. Don’t think a second that I don’t have massive respect for @certik and others who are doing something about fortran instead of the others just endlessly praising a de-facto dead language and doing nothing about it. However, that also doesn’t mean that I agree with his outlook on Fortran. I think the current team is exceedingly productive and the ecosystem is ripe for a fast Rust-like array language. So in my opinion if they drop off all the, excuse my french, historical shit behind and start tabula rasa, they would be adopted much better. And clearly they have the knowledge. I don’t quite understand the loyalty to this lang. This is not to belittle anything, (other than the users I mentioned above), people are even trying to leave C/C++ behind. So if you think fortran is worth keeping imagine C++. Hence nothing is sacred and if and when the times come, we should be objective. Not expecting a committee to come and save us.
Folks here didn’t quite get the point but let me repeat, first of all matlab scripting language uses 1 indexing because of legacy but internally uses 0 indexing so remove that one off the list. Then you have fortran, julia and R which are practically followers of fortran. So it is fortran’s fault to start with. If you are still at this debate, lots of surprises waiting for you when you want to link to other languages and be adopted for the low level programming.
We touched this issue a few times over scipy when there was some activity about @zaikunzhang 's PRIMA rewrite. Judging by the avatars, I think it was @hkvzjal who were looking into the example I gave but I don’t know what happened afterwards about noncontiguous array slices of NumPy. If you keep using contiguous column major arrays, you can kiss goodbye to your precious performance gains if you glue to anything that is row major. Since fortran and julia are the ones left in the column major club, you would be playing on your own. If you haven’t done so, work a bit with C code or NumPy or Rust and you’ll see that more often than not you would be creating array copies left and right because of the mismatches. So I’m not making this up from my ass. That’s all we do in the Cythonized code dealing with this silly transposing interplay.
Second to last, because I can’t miss that, Pearu Peterson et al. did not adopt F77. He is the author of f2py
that allowed fortran to be relevant for all python users for decades. So I’d careful who I am pointing to as not being up-to-par fortranners. They included these F77 libraries because, there was no other better version, in fact there was nothing else. Like I said we should not have been doing these rewrites. People like you who give a shit about fortran should have been doing this out of historical responsibility. Also you can feel the fortran/matlab influence in NumPy. So I think you should thank them the most for some of the choices made in NumPy even though there is nothing in NumPy that is fortran (including the blas/lapack). So basically they are doing this as a favor to the community.
@ivanpribec writes
Unfortunately, I cannot live off of refactoring Fortran codes, hence I had no motivation to contribute it to SciPy, despite being a user
None of us are living off of this. I know you mean well, I truly do. But that’s the thing I see here in this forum in the last few days. Because I read a lot in the last few days We are treated as if we are some sort of exclusive money making snob asshole club. One guy was angry because I used “get rid of” in a sentence and wrote off the whole SciPy as shit, I mean that was fun Some guy got triggered because I said incompetent although the definition fits perfectly like a glove (there is a thread about an http library next to this one I mean all the glory but did you folks write any nonscientific code lately? Do you realise how much fortran is missing?). Someone else thinks Julia people are trying to overthrow Fortran like as if it is not self-destructing itself. Apologies for saying this but this is just hilarious. Never had this much interesting takes in a forum for a while.
Someone was arguing that HPC is the only thing that matters in scientific code and LANL people are just delusional just to question Fortran usage. The audacity! The blasphemy! Some think we are bunch of amateurs in other languages punching shitty code and having wrong ideas about how scientific algorithms should work. Whether that is true or not still, I have to confess, I find it really funny but also more importantly disheartening. Some of the maintainers of these Python packages are actually experts on their fields. They are all sacrificing from their social life to do a little bit of thankless good only to get shit like above . Don’t make the noob mistake that just because we don’t like fortran, we don’t know how to work with it. Like I said I started with matlab 5 and fortran 90. I’m not a good fortran coder, and I would like to stay that way. Incredibly enough I never thought that keywords can be variables. Because which competent language allow that?
I hope this gives some insight to the brain of a single maintainer of an obscure package in python who does not represent that package or python or its users. Don’t take things this seriously. Only mistake I did was to give a candid response about fortran because I do think it is just unsalvagable mess in 2023. But hey it is my candid stupid opinion. It is you that has the option to get offended
Anyways I spent enough time here, though it was fun has to come to an end. I thank again all for their explanation and help. Take care. Looking forward to what comes out from LFortran folks.
P.S.: Julia included, many scientific packages are rewriting their own {GE/SY/PO}TRF
because small array performance n ~< 100-150
of LAPACK libs MKL/OpenBLAS/BLIS is not good enough and have large overhead (nobody uses reference implementation from netlib anyways which is supposed to be fast because fortran right?). So I would like to emphasize, nothing is sacred or untouchable.