SciPy: Interested in adopting PRIMA, but little appetite for more Fortran code

https://fortran-lang.discourse.group/t/scipy-there-is-interest-in-adopting-this-prima-implementation-but-little-appetite-for-taking-on-more-fortran-code/5821

by zaikunzhang • 3 years ago

126 68 3 years ago

fortran-lang.discourse.group

zaikunzhang10• 3 years ago

see also

Optimization Without Derivatives: PRIMA Fortran Version and Inclusion in SciPy, https://news.ycombinator.com/item?id=35959991

SciPy enhancement: The Fortran 77 implementation of COBYLA is buggy and challenging to maintain. Switch to the PRIMA implementation? #18118, https://github.com/scipy/scipy/issues/18118#issuecomment-155...

bouchard3• 3 years ago

I understand the lack of appetite for old FORTRAN 77 code, but PRIMA is a modern Fortran implementation.

You only need Python's ctypes module to load the compiled Fortran library thanks to Fortran's iso_c_binding module.

Using f2c as suggested seems like a worse alternative, not a solution.

wirrbel• 3 years ago

Using f2c as suggested seems like a worse alternative, not a solution.

essentially madness.

Modern fortran is quite nice, much nicer than C from a Python programmer's perspective.

zaikunzhang10• 3 years ago

Modern fortran is quite nice, much nicer than C from a Python programmer's perspective.

Totally agree. I hope everyone making comments on Fortran can spend a few minitues reading what modern Fortran is, e.g., by checking

https://fortran-lang.org/en/

drbaba• 3 years ago

I used Fortran 2008 a lot some years ago; initially I wanted to write a C++ wrapper around a Fortran library that seemed useful, but ended up liking the language better than C++ and just did the whole project in it.

I think a lot of people don’t realize that modern Fortran actually supports object-oriented programming and Matlab-style array slicing, and comes with built-in support for e.g. matrices and complex numbers. It even has some features inspired by functional languages, like being able to mark functions as pure (side-effect free). Combined with its general speed, this makes it a pretty nice language for numerical code.

Although for most new users, I think Julia might be a better option than Fortran these days.

noobermin• 3 years ago

This reads like the maintainer, mind you a maintainer of scipy, the premier foundational library in numerical python, does not even understand what PRIMA is or what modern fortran is like. I get people are busy, but this shakes my opinion of scipy a little bit. It doesn't take that long to read the readme... It does not reflect well on them.

fluidcruft• 3 years ago

I've been trying to debug some of scipy's statistics and it's... incomprehensible. I am not impressed. I'll happily take fortran over its growing morass of wierdo boost wrappers.

tomr_stargazer• 3 years ago

I thought that the PRIMA developer’s comment in this thread was really mature and thoughtful. Reproduced below (minus some formatting):

My goal is to make Professor Powell 12’s solvers as accessible as possible to scientists, engineers, and algorithm researchers. I am not particularly in favor of or against any language. I hope that everyone can easily use Powell’s solvers in her/his favorite languages.
The first implementation of PRIMA is in modern Fortran simply because Powell’s implementation was in Fortran 77. Using Fortran, I can systematically verify the bit-to-bit faithfulness of the modernized implementation (not only “faithful up to an epsilon”). In addition, the intrinsic support for matrix-vector calculations is a strong advantage when developing reference implementations (or templates) of numerical solvers — most numerical algorithms are combinations of such calculations anyway.
The major motivation for developing the modern Fortran version is to provide a reference for the implementation in other languages, namely Python, MATLAB, C++, Julia, and R. A reference implementation must be structured, modularized, readable, understandable, and extendable. The original Fortran 77 code is a true masterpiece, but it is not proper at all for being used as a reference implementation. You do not want to use a spaghetti-style codebase with 244 GOTOs as a reference, or your implementation will be of the same style.
Putting it more straightforwardly, I implemented the modern-Fortran version of PRIMA in order to develop versions that are entirely Fortran-free 3.
Coming back to the point, I fully understand why the SciPy community has “little appetite for taking on more Fortran code”. The reputation of Fortran has been damaged over the years. I do not agree with the damaged reputation, and I feel sorry for those who do not have a chance to know (or refuse to know) modern Fortran due to this false reputation, but I do not blame them. It is the responsibility of the Fortran community to re-establish the reputation. I do not regard it as a bully to request non-Fortran implementations. It is not a question of surrendering or not.
I do consider myself a member of the Fortran community. Taking my share of the aforementioned responsibilty, I will try to promote the usage of modern Fortran via the PRIMA project. There is nothing more convincing than a successful real-life project.
For the inclusion of PRIMA in SciPy, I will keep communicating with both the SciPy and the Fortran community (e.g., those on this discourse), trying to find the best route. As pointed out by others, f2c is not an option due to its incapability of handling modern Fortran. Official C++ and Python implementations are being planned, but they will not be delivered in the near future. The most probable and practical solution, as suggested here and under the other thread, is to wrap the modern Fortran implementation of PRIMA using iso_c_bindings + ctypes or similar facilities. I hope the SciPy maintainers will accept this solution.

counters2• 3 years ago

Yeah, I agree 100% - I really support the approach of the developer here and am totally aligned with all the reasoning.

The libprima/prima codebase is very readable, even if you're not accustomed to modern Fortran (let alone Fortran period). It has fantastic comments throughout the numerical algorithms, and even though there are a lot of lines of code, in most places it really seems to be a minimally-complex implementation with very little magic. I haven't built anything against the codebase yet but based on the examples I feel it will be far easier than many other libraries out there.

A motivated developer could very quickly port this code to their preferred numerical programming language.

zaikunzhang10• 3 years ago

Thank you so much for your very warm and encouraging comments!

PRIMA has been a black hole that absorbs all my time and energy in the past three years, which even puts my career (as a junor professor) in real danger. The positive feedback like yours is vital for me. Without it, I would not have the energe or courage to continue. Many thanks!

A motivated developer could very quickly port this code to their preferred numerical programming language.

I am much glad to hear a person other than myself saying this. It is the very reason why I develop this reference implementation. PRIMA achieves its success if others can implement Powell's solvers to high quality using PRIMA as a referene, without the genius like Powell and without the long-term experience and strugling like me.

Many thanks! --- Zaikun

zaikunzhang10• 3 years ago

Thank you very much for quoting this! --- Zaikun

radarsat1• 3 years ago

I wonder why not just try to use more Pythonic solutions such as numba within scipy. I realize that maybe not in all cases can the same performance as Fortran be achieved but it strikes me as much more hackable for Python programmers than having to learn an new (old) language just for certain core algorithms. One of the principal reasons for switching to C implementations is because more people understand C than Fortran. This seems doubly true for Python-based jit solutions like numba or jax, so why not just skip C altogether and go for a more native port?

adastra226• 3 years ago

Numpy is Fortran btw.

remram• 3 years ago

From numpy.org:

The core of NumPy is well-optimized C code.

whimsicalism8• 3 years ago

no?

bee_rider6• 3 years ago

Numpy is most importantly BLAS and LAPACK, which can be implemented in any language. The best CPU version is probably still MKL, which is a mix of Fortran, C, and assembly. C is used for these kinds of projects because it is easy to sprinkle assembly into C.

Fortran solves a particular problem: produce a pretty good binary from readable code that a normal PHD student or engineer can write. This is not really the problem that a tuned BLAS/LAPACK library has to face, vendors throw very experienced folks at these libraries to chase after those last 5% performance boosts that Fortran leaves on the table.

zaikunzhang10• 3 years ago

A native port is indeed planned. However, since we are talking about a project of about 10K lines of code (though many of the lines are comments), such a port will not be delivered very soon.

In fact, native implementations of PRIMA in Python, MATLAB, C++, Julia, and R will all be done in the future. See https://github.com/libprima/prima#other-languages . But they will take time. PRIMA has been a one-man project since it started three years ago. Community help is greatly needed.

Thanks.

adastra226• 3 years ago

This is crazy. Get over yourself! Whether you like it or not, Fortran is still THE industry standard for scientific codes, and the fastest compiled language on many architectures.

raoulj• 3 years ago

In the defense of the SciPy response, while there may be some programming language zealotry behind the position, it seems to be from a lack of expertise in SciPy to maintain its Fortran code that's pushing SciPy to avoid new Fortran additions. Hopefully OP gets PRIMA through! Been admiring this effort, Dr. Zhang.

ok1234562• 3 years ago

F2py, which is already a standard part of numpy, does the heaving lifting and generates idiomatic python wrappers for Fortan code.

It's just zealotry of them going "eeew fortran that's like really really old". Fortran got multidimensional dense numeric arrays right. Modern Fortran allows you to do elemental/vectorized operations on them. C does not even come close, and requires tons of preprocessor macros if you try. Every other 'modern' language introduces a bunch of tooling and dependencies just to get where Fortran is, because they just don't have a focus on scientific computing.

zaikunzhang10• 3 years ago

Thank you so much for your encouragement! Such positive feedback means a lot to me.

cozzyd3• 3 years ago

As long as PRIMA is maintained (which seems likely independently), it's only the bindings that need to be maintained. And it's trivial to write a Fortran interface in C or C++ (perhaps that can be included in PRIMA anyway) and then SciPy can pretend like it's not Fortran...

asdfman123• 3 years ago

Are you saying they're acting like PRIMA-donnas?

travisporter• 3 years ago

Who is this insane developer who links twitter to every text selection on this page!@

cozzyd3• 3 years ago

Looking at the original issue, scipy maintainers don't want to use LGPL code (why?!?) but they use code that appears in Numerical Recipes (which has super onerous licensing)?!?

jabl• 3 years ago

First, if you read through that scipy issue (https://github.com/scipy/scipy/issues/18118 ) the author was willing and able to relicense PRIMA under a 3-clause BSD license which is perfectly acceptable for scipy.

For the numerical recipes reference, there is a mention that scipy uses a slightly improved version of Powell's algorithm that is originally due to Forman Acton and presumably published in his popular book on numerical analysis, and that also happens to be described & included in numerical recipes. That is, unless the code scipy uses is copied from numerical recipes, which I presume it isn't, NR having the same algorithm doesn't mean that every other independent implementation of that algorithm falls under NR copyright.

cozzyd3• 3 years ago

thanks, I misunderstood the NR reference. (In my field... it is common for people to copy and paste terrible NR code without understanding how bad it is or the licensing terms...).

domcoltd• 3 years ago

Hopefully, the SciPy community can stay open-minded about modern Fortran libraries.

Modern Fortran is quite different from Fortran 77, while being as powerful, if not more.

In addition, there has been a significant community effort on improving and modernising the legacy packages, the ecosystem, and the language itself.

With projects like LFortran (https://lfortran.org/), fpm (https://github.com/fortran-lang/fpm), and stdlib (https://github.com/fortran-lang/stdlib), I believe that Fortran will enjoy prosperity again.

zaikunzhang10• 3 years ago

Each of LFortran, fpm, and stdlib deserve a HN front-page story :)

nvy• 3 years ago

Any resources for learning Fortran in 202X that don't bother going over the old cruft?

Last time I looked (admittedly, years ago), most resources seemed to be in the vein of "here's Fortran 1990, and once you've learned that here are a few chapters on F03 and F08, treated separately".

It'd be great to learn there exists a "Fortran: the good parts" sort of thing that just completely elides the outdated stuff.

bouchard3• 3 years ago

I found "Modern Fortran: Building Efficient Parallel Applications" to be a pretty good read. See https://milancurcic.com/book/

Also https://fortran-lang.org/en/ is a good hub to find more information, packages, etc.

nunuvit2• 3 years ago

https://fortran-lang.org/en/learn/quickstart/

zaikunzhang10• 3 years ago

A list of books, tutorials, and online courses:

https://fortran-lang.org/en/learn/

BTW, the Fortran Lang website (https://fortran-lang.org/) is THE official website of modern Fortran, maintained by the modern Fortran community.

You may also take a look at the modern Fortran discourse at

https://fortran-lang.discourse.group/

where you can find all the gurus of modern Fortran. They are very willing to help if you post a question. I have asked all my Fortran questions over there during the development of PRIMA.

whinvik• 3 years ago

Even though it's not great but I understand why this maybe the direction SciPy maintainers want to go in. They have limited capacity and I think focusing on a single language makes sense especially if there are way more developers in that language, C.

One way to solve this would be of course, to have more people from the Fortran community being part of the SciPy community. But I don't know how big the Fortran community really is to be able to do that.

nunuvit2• 3 years ago

There just isn't much interest. HPC focuses on distributed and product doesn't care about open source. That mostly leaves just professors and oddballs.

ChrisRackauckas• 3 years ago

Interesting response. I develop the Julia SciML organization https://sciml.ai/ and we'd be more than happy to work with you to get wrappers for PRIMA into Optimization.jl's general interface (https://docs.sciml.ai/Optimization/stable/). Please get in touch and we can figure out how to set this all up. I personally would be curious to try this out and do some benchmarks against nlopt methods.

SciPy: Interested in adopting PRIMA, but little appetite for more Fortran code

Related Stories